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CHAPTER  1 


INTRODUCTION 


Fundamental  to  positioning  is  the  concept  of  height  or  elevation  above  some  vertical 
reference  datum.  Since  the  process  of  determining  heights  (or  height  differences)  using  spirit 
leveling  demands  that  the  vertical  reference  surface  be  an  equipotential  surface,  historic^ly  the 
geoid  has  been  used  for  ^at  purpose.  By  its  very  definition,  the  geoid  refers  to  an  equipotential 
surface  of  the  gravity  field  of  the  earth  that  most  closely  approximates  the  mean  sea  level  in  the 
spatial  least  squares  sense,  at  a  given  time. 

There  are  more  than  three  hundred  regional  vertical  datums  in  the  world  today,  most  of 
them  are  based  on  the  long  term  mean  value  of  the  local  mean  sea  level  observed  at  one  or  more 
primary  tide  gauge  stations,  with  an  implicit  assumption  that  the  mean  sea  level  at  these  tide 
gauge  sites  coincide  with  the  geoid.  With  the  present  day  accuracy  requirement  of  ±10  cm  for 
the  vertical  datum  definition  (Rapp,  1983),  the  mean  value  of  the  local  mean  sea  level  observed 
at  the  tide  gauge  sites  cannot  be  consider^  to  coincide  with  the  geoid.  The  departure  of  the 
mean  sea  level  (with  respect  to  a  specific  time  period)  from  the  geoid  is  termed  as  Sea  Surface 
Topography  (SST),  which  has  a  magnitude  of  up  to  ±2m  globally  even  over  long  time  scales. 
SST  can  be  considered  as  a  function  of  the  geographic  location  (<|),  X)  and  the  time  (t).  Therefore 
it  varies  temporally  and  spatially.  Most  of  ^e  temporal  variations  within  a  specific  period  can  be 
averaged  out  using  appropriate  filters,  but  those  which  are  not  cyclic  in  nature  must  be 
consigned  for  the  computation  of  SST  at  the  epoch  for  which  the  mean  sea  level  is  cmnputed. 
Hiis  leads  to  the  concept  of  an  ep<xh-dependent  datum  (Cross  et  al.,  1987).  This  datum  is 
established  by  taking  tide  gauge  readings  at  one  or  more  tide  ^uge  stations  over  a  particular  time 
span  and  reducing  these  measurements  to  a  specific  epoch  within  the  time  interval.  All  heights 
are  then  referred  to  this  adopted  vertical  reference  datum.  Mean  sea  level  at  a  later  epoch  will 
not  coincide  with  this  reference  surface. 

Since  the  vertical  datums  adopted  by  different  countries  are  based  on  mean  sea  levels 
computed  for  different  epochs,  using  different  procedures,  and  tide  gauge  measurements,  they 
^ffer  fiom  each  other  as  well  as  wi^  the  geoid  to  an  extent  of  ±2m.  Estimated  vertical  datum 
inconsistencies  5h  (Laskowski,  1983)  and  its  propagated  effects  on  terrestrial  gravity  anomalies 
over  continental  plates  are  given  in  Heck  (1990,  Table  1).  The  magnitude  of  the  effects  of 
vertical  datum  inconsistencies  on  gravity  antnnalies  are  on  the  order  of  ^.5  mgal  (ibid.).  We  see 
this  effect  reflected  as  systematic  regional  biases  in  gravity  anomaly  data  banks. 

A  few  years  back,  geodetic  applications  were  not  affected  by  such  discrepancies  in 
vertical  datums  as  they  have  played  a  small  role  in  the  problems  of  geodesy.  As  greater  accuracy 
is  required  and  more  accurate  measurements  are  made  now,  the  problems  caused  by  vertictU 
datum  inconsistencies  have  started  surfacing  more  dramatically  than  before  (Rapp,  1987).  The 
requirements  for  defining  and  realizing  a  global  vertical  datum  with  a  precision  on  the  dm-  or 
even  cm-level  has  become  essential  now,  nxne  than  before,  due  to  the  improved  measurement 
accuracy  provided  by  the  space  geodetic  techniques.  Several  arguments  emphasizing  the  need 
for  a  global  vertical  datum  have  been  discussed  by  Heck  and  Rummel  (1990)  and  Carter  et  al. 
(1989).  A  few  of  the  strong  arguments  are  as  follows; 
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•  The  establishment  of  a  global  vertical  datum  will  provide  a  means  to  accurately  connect 
national  and  continental  vertical  datums. 

•  Removal  of  the  systematic  regional  biases  in  gravity  anomaly  data  banks,  noted  earlier  due  to 
referencing  heights  to  different  level  surfaces  in  different  vertical  datum  zones,  requires  the 
definition  of  a  global  vertical  datum,  so  that  the  gravity  anomalies  can  now  refer  to  one 
unique  geopotential  surface. 

•  Comparisons  between  the  results  of  geodetic  leveling  and  oceanographic  procedures  for 
determining  sea  surface  slopes  over  long  distances  require  a  consistent  vertical  datum. 

In  addition  to  the  above  arguments,  the  rapid  growth  in  positioning  through  the  Global 
Positioning  System  (GPS)  has  created  a  need  for  better  understanding  of  orthometric  heights 
globally  and  how  they  can  be  determined  using  GPS. 

1.1  Eroblcm  and  the  SackgTQund 

Attempts  to  come  up  with  a  global  vertical  datum  or  a  world  vertical  network,  either  by 
connecting  the  existing  regional  vertical  networks  or  by  defining  a  global  datum  irrespective  of 
any  speciHc  regional  vertical  datum,  is  not  new  to  the  geodetic  community.  A  special  study 
group  1.75  on  the  World  Vertical  Reference  System  was  established  at  the  1983  Hamburg 
lUGG/IAG  meeting,  with  the  primary  purpose  to  study  methods  for  the  determination  of  a 
unique  vertical  reference  surface  that  could  be  used  on  a  global  basis,  under  the  chairmanship  of 
Prof.  Richard  H.  Rapp  of  The  Ohio  State  University.  T^e  study  group  recommended  that  the 
need  to  study  the  vertical  datum  problem  is  more  relevant  in  this  space  geodetic  age  and  we 
should  take  advantage  of  the  new  data  and  programs  that  are  now  on  the  horizon  (Rapp,  1987). 
“The  First  International  Conference  on  Geodetic  Aspects  of  the  Law  of  the  Sea”  organized  by  the 
lAG  Subcommission  GALOS,  at  Denpasar,  Bali-Indonesia  8-13  June  1992,  by  their  resolution 
no.  3  noted  the  imponance  of  consistency  in  defining  sea  coast  positions  for  boundary  puiposes, 
urged  the  International  Association  of  Geodesy  to  pursue  with  other  interested  organizations  the 
possibility  of  advocating  a  global  vertical  datum  and  also  to  investigate  a  suitable  definition  of 
such  a  datum. 

In  addition  to  several  International  commission  reports  and  seminar  proceedings, 
numerous  papers  have  also  been  published  in  the  past  twenty  years  to  document  the  problem  of 
the  global  vertical  datum  definition  and  proposed  solutions,  llie  relative  merits  and  ^merits  of 
the  oceanographic  approach  suggested  by  Cartwright  (1985),  Altimetric-Gravimetric  approach 
by  Mather  et  al.  (1978),  and  Integrated  Geodesy  approach  of  Hein  and  Eissfeller  (1985)  are 
discussed  in  detail  by  Pavlis  (1991)  and  partially  by  Heck  and  Rumroel  (1990). 

With  respect  to  the  global  vertical  datum  definition,  the  combination  of  geometrical  and 
geodetic  data  is  proposed  in  many  ways.  Colombo  (1980),  Ra^J  (1980),  Hajela  (1985)  and 
Pavlis  (1991)  have  suggested  techniques  for  venical  datum  connection  by  estimating  the 
geopotential  differences  over  intercontinental  locations  using  satellite  and  terrestrial 
measurements.  They  use  three-dimensional  geocentric  coordinates  of  at  least  one  fundamental 
point  in  each  vertical  datum  system.  The  geoid  height  at  these  locations  can  be  calculated  from 
some  high-degree  geopotentitd  model,  and  the  influence  of  local  differences  between  the  actual 
gravity  freld  and  the  global  model  is  determined  from  gravity  anomalies  given  in  some  spherical 
cap  around  the  fundamental  point.  The  estimated  accuracy  of  vertical  datum  connection  between 
different  continents  is  about  ±15  to  ±17  cm  (Pavlis,  1991). 

Another  approach  described  by  Rummel  and  Teunissen  (1988),  Heck  and  Rummel 
(1990),  and  Xu  and  Rummel  (1991),  uses  globally  distributed  gravity  anomalies.  These  values 
are  given  referring  to  different  vertical  datum  realizations  for  different  parts  of  the  eanh.  Fixing 
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one  datum  as  reference,  all  the  others  may  be  determined  if  the  onhometric  heights  of  the 
stations  in  local  vertical  datum  together  with  their  geocentric  coordinates  are  known  for  at  least 
one  fundamental  benchmaiic  in  each  vertical  datum  zone. 

The  ideas  described  in  the  above  two  paragraphs  clearly  show  that  the  studies  so  far  made 
lead  to  either  a  connection  between  selected  vertical  datums  or  suggest  a  procedure  to  determine 
the  potential  difference  between  selects!  benchmarks  on  two  different  regional  datums  which 
can  lead  to  a  vertical  datum  connection.  The  concept  of  defining  and  realizing  a  global  vertical 
datum,  independent  of  any  particular  regional  datum,  has  been  discussed  by  Rapp  and 
Balasubramania  (1992).  This  dissertation  woric  enhances  the  ideas  formulated  in  the  report 
(ibid.)  to  come  up  with  a  practically  realizable  global  vertical  datum. 

1.2  Organization  and  Scope  of  the  Present  Study 

From  the  above  sections  we  see  that  most  of  the  countries  in  attempting  to  define  a 
vertical  reference  surface  have  relied  on  local  mean  sea  level  as  the  surface  that  is  accessible 
through  tide  gauge  measurements.  Since  we  are  now  moving  into  an  era  where  more  precise  and 
more  consistent  height  datum  determinations  are  essential,  there  is  a  necessity  to  come  up  with  a 
vertical  datum  that  can  be  used  on  a  global  basis.  There  are  basically  two  issues  involved  in  this 
global  vertical  datum  development,  that  is,  connection  of  existing  vertical  datums  and  deHnition 
of  a  global  datum  independent  of  any  particular  regional  datum.  As  mentioned  earlier,  these 
issues  have  been  discussed  in  the  literature  in  the  past  few  years,  but  no  specific  action  has  been 
taken  so  far. 

This  dissertation  will  attempt  to  address  these  issues.  In  addition  to  detailed 
mathematical  modeling  to  define  and  realize  a  global  vertical  datum  to  an  accuracy  of  ±10cm 
after  considering  all  the  error  sources  in  the  observed  data,  a  test  calculation  has  also  been 
included  to  come  up  with  the  first  iteration  global  vertical  datum  to  an  accuracy  permitted  by  the 
available  data  today. 

Two  different  approaches  are  developed  fOT  this  purpose.  The  first  approach  uses  the  best 
available  station  positions  from  the  jmsent  and  future  space  geodetic  networks,  highly  accurate 
geopotential  models  and  surface  gravity  data  around  the  space  stations  for  defining  and  realizing 
a  global  vertical  datum.  This  approach  can  be  considered  as  an  ideal  solution  of  the  problem,  but 
needs  to  be  treated  as  a  long  term  proposal  since  the  data  required  are  not  available  all  over  the 
world  and/or  to  required  precision. 

The  other  approach  can  be  considered  more  of  an  operational  development  for  immediate 
implementation  but  less  accurate  than  the  first  approach.  This  approach  is  mainly  dependent  on 
the  widely  available  space  geodetic  data  from  (3PS  (Global  Positioning  System)  and  DORIS 
(Doppler  Orbitography  and  Radiopositioning  Integrated  by  Satellite)  tracking  networks  and 
highly  accurate  regional  geoid  height  models. 

The  dissertation  structure  is  as  follows; 

Chapter  2  explains  the  height  datum  fundamentals,  bringing  out  the  different  types  of 
heights  that  are  in  use  today  with  their  definition,  computation  and  intercomparisons.  Procedures 
for  establishing  a  more  accurate  and  meaningful  regional  vertical  datums,  which  in  turn  can  help 
in  realizing  an  accurate  global  vertical  datum,  are  also  discussed  in  this  chapter. 

The  primary  step  in  the  development  of  global  vertical  datum  involves  generating 
mathematical  models  and  defining  required  observations.  Since  the  data  availability  in  the 
various  parts  of  the  world  is  not  uniform,  the  development  of  mathematical  models  under 
different  data  scenarios  is  discussed  in  Chapter  3.  Consideration  of  different  error  sources  in  the 
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available  terrestrial  gravity  anomalies  and  their  impact  on  computed  gravimetric  undulations  are 
included  in  this  chapter.  Also  included  are  the  observation  equation  model  development  and  the 
adjustment  procedure  to  estimate  the  parameters  defining  the  global  vertical  datum. 

In  Chapter  4,  a  test  calculation  for  implementing  the  Global  Vertical  Datum  based  on  the 
available  data  on  six  of  the  regional  vertical  datums  is  included.  The  six  regional  vertical  datums 
considered  in  the  test  calculations  are  viz..  North  American  Vertical  Datum  (NAVD  88), 
Odnance  Datum  Newlyn  (ODN,  UK),  IGN  69  prance),  NN  (Germany),  AHD:71  (Australia) 
and  Scandinavian  datums.  Though  there  is  no  single  regional  vertical  datum  as  Scandinavian 
datum,  the  individual  datums  RH70  (Sw^n),  N60  (Finland)  and  NN1954,  NNN1957  (Norway) 
were  combined  to  serve  as  a  single  regional  datum  for  computational  convenience.  More  details 
on  the  Scandinavian  height  datum  are  given  in  Chapter  4.  Two  different  computational 
procedures,  namely.  Least  Squares  Collocation  technique  and  modifled  Stokes’  integration 
technique,  have  been  used  to  compute  the  gravimetric  undulation.  These  results  have  been 
compa^  with  undulation  values  from  regional  geoid  height  models  (for  U.S.  and  Scandinavian 
space  geodetic  stations  only)  which  are  generally  developed  using  Fast-Fourier  Transform  (FFT) 
techniques. 

Chapter  S  deals  mainly  with  the  simplified  approach  in  realizing  a  global  vertical  datum. 
The  theory  behind  the  model  development  along  with  test  results  and  plots  in  the  area  of  GPS 
traverses  in  the  U.S.A.,  Canada  and  Western  Europe  are  included  in  this  chapter.  Height  bias 
modeling  has  been  done  in  four  test  areas  in  Canada  and  the  U.S.A. 

While  Chapter  6  brings  out  the  summary  and  conclusions  drawn  from  this  dissertation 
work.  Chapter  7  lists  the  recommendations  for  the  acquisition  and  use  of  data  for  defining  and 
realizing  a  global  vertical  datum  to  a  precision  of  a  few  centimeters. 
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CHAPTER2 


PRELIMINARIES  AND  DEFINITIONS 


There  are  fundamentally  two  different  height  systenos  available  today  which  can  be 
realized  with  approximately  the  same  order  of  accuracy  (Torge,  1987).  The  first  one  is  the 
gravity  field  related  hei^t  system  which  is  based  on  leveling  and  gravity  data  alcmg  the  leveling 
lines.  These  observations  lead  to  get^tential  heights,  if  referred  to  a  zero  height  level  surface. 
The  other  one  is  the  ellipsoidal  height  system  realized  through  satellite  techniques  such  as  SLR, 
GPS-positioning  etc.  Original  results  are  the  cartesian  coordinates  of  the  stations  in  a  near 
geocentric  coordinate  system  which  can  be  transfnred  to  heights  called  ellipsoidal  heights  over 
the  reference  level  ellipsoid  Knowing  the  gravimetrically  derived  geoid  information  from  geoid 
height  nxxlels,  the  orthometric  height  of  the  station  can  be  computed  in  very  good  approximation 
using  the  relation: 


H  =  h-N  (2.1) 

In  the  above  equation,  the  ellipsoidal  height,  h,  is  the  distance  from  a  defined  ellipsoid  to  a  point, 
measured  along  the  straight  line  that  is  peipendicular  to  the  ellipsoid  and  passes  through  the 
point.  The  orthometric  height,  H,  is  the  separation  between  the  geoid  and  the  point  measured 
along  the  curved  vertical  defrn^  by  the  equipotential  surfaces  of  dte  earth’s  gravity  field  and  N, 
the  geoid  undulation,  is  the  separation  between  the  defiited  ellipsoid  and  geoid  (Rapp,  1992). 

The  relative  orthometric  height  accuracy  that  can  be  achieved  using  the  above  two  height 
systems  are  given  by  Niemeier  (1987,  Table  3)  and  shown  in  Table  2.1. 

Table  2. 1  Height  System  and  Accuracy  For  Cottimmily  Used  Observation  Techniques 


Height 

system 

Observation 

technique 

Additional 

information 

Quantities  available 
for  ortho,  ht 
computation 

Relative  accuracy 

of  computed 
orthometric  hts. 

iSSiMSHi 

Surface  gravity 

height  increments 
potential  differences 

IIP9IIII 

tteight  increnaents 
potential  differences 

geoid 

ellipsoidal  height 
differences 
geoidal  undulation 
^fferences 

Though  several  papers  have  been  published  on  the  topic  of  obtaining  the  mthometric 
heights  using  GPS  and  other  space  geodetic  techniques  in  the  recent  years,  no  county  in  the 
world  has  so  far  ad(q)ted  the  t^ometric  heights  obtained  from  space  geodetic  techniques  for 
defining  their  vertical  network.  The  main  reason  for  this  is  the  insufficient  knowledge  of  the 
regional  geoid  heists  to  the  accuracy  of  leveled  heights.  With  constant  efforts  to  improve  the 
regional  geoid  height  models  and  also  the  ellipsoidal  height  accuracy  for  GPS  observatitHis,  it 
may  be  possible  in  the  future  to  replace  conventional  leveling  by  satellite  based  techniques.  But 
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as  of  today,  since  the  conventional  leveling  is  the  only  technique  accepted  by  the  geodetic 
community  for  defining  precise  heights,  the  following  sections  discuss  the  various  heights  and 
height  systems  that  are  used  in  different  parts  of  the  world  today  which  are  based  mainly  on  the 
‘gravity  related  height  system’. 

Detailed  discussions  on  height  systems,  their  measuring  prowdures  and  the  treatment  of 
height  networks  may  be  found  in  the  books  of  Heiskanen  and  Moritz  (1967,  Chapter  4),  Torge 
(1991,  Chapter  6),  Vanicek  and  Krakiwsky  (1986,  Chapter  19),  Bomfo^  (1980,  Chapter  3)  and 
others.  Additional  details  on  various  kinds  of  heights  and  their  computation  may  be  found  in 
Krakiwsky  (1965)  and  exhaustive  treatise  only  on  onhometric  heights  in  Rapp  (1961). 

2d  Height  Definitiops  and  Height  Systems 

For  decades  the  well  known  geometric  or  spirit  leveling  was  the  only  technique  used  for 
determination  of  precise  heights.  New  leveling  techniques  such  as  nmtorized  geometric  leveling, 
motorized  trigonometric  leveling  and  motod^  3D-traversing,  etc.,  have  been  developed  now, 
to  achieve  the  same  accuracy  as  conventional  leveling  but  to  reduce  the  time  and  expenses 
involved  to  perform  convendonal  leveling  for  extended  networks  (Niemeier,  1986).  All  these 
measurements  (leveling,  theodolite  measurements,  etc.)  are  almost  exclusively  referred  to  the 
system  of  level  surfaces  and  plumb  lines,  the  geoid  playing  an  essential  part. 

2.1.1  Level  Surfaces  and  Geopotential  Numbers 

Level  instruments  and  theodolites,  by  virtue  of  their  dependence  on  the  local  horizon  with 
respect  to  their  orientation  by  the  actual  gravity  vector  "f,  are  influenced  by  the  actual  earth’s 
gravity  field.  The  surfaces 

W(x,  y,  z)  =  Wo  =  const.,  (2.2) 

on  which  the  potential  W  is  constant,  are  called  equipotendal  surfaces  or  level  surfaces.  Due  to 
the  figure  of  the  earth,  the  earth  rotation  and  inhomogeneous  mass  distribution  the  level  surfaces 
near  the  earth’s  surface  are  not  parallel  and  the  lines  that  intersect  all  equipotential  surfaces 
orthogonally  are  not  exactly  straight  but  slightly  curved  (refer  Figure  2.1). 
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For  a  displacement  along  the  level  surface,  the  potential  difference  dW  »  0,  for  displacement 
along  die  outer  surface  ncnrmal  n  the  potential  difmence  is 


dW  =  -gdn  (2.3) 

where  dn  is  the  separation  between  the  level  surfaces  and  g  the  magnitude  of  gravity  at  the 
plumbline  (see  Figure  2.3  also),  ^uadon  (2.3)  provides  the  link  between  the  potential 
difference  dW  and  Ae  difference  in  height  of  neighborag  level  surfaces. 

As  reference  for  most  height  systems  the  gecnd  plays  an  inqiortant  role.  Connected  with 
the  geoid  is  the  zero  p<^t  or  datum  problem.  C3«e  has  to  define  which  point  Pq  is  located  on 
geoid  having  a  potential  Wq.  Since  it  was  tacidy  believed  that  the  mean  sea  level  should 
Aeoretically  coincide  with  the  geoid,  or  that  the  departure  of  die  two  surfaces  was  negligible  (or 
unknown),  locating  the  geoid  or  fixing  the  point  Po  with  respect  to  a  reference  benchmark  on  Ae 
shore  was  reduced  to  the  task  of  determining  the  position  of  die  mean  sea  level.  In  most  of  the 
regional  vertical  datums  in  the  world  today,  the  point  Pq  has  been  identified  to  lie  on  the  long 
term  mean  value  of  the  local  mean  sea  level  (generally  with  18.6  years  of  sea  level  data)  with 
assi^ed  potendal  W  =  Wq.  Establishment  of  a  reference  benchmark  with  respect  to  Po  is  shown 
in  Figure  2.2  (given  as  Figure  19.1,  Vanicek  and  Kraldwsky,  1986). 
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Figure  2.2  Establishment  of  Height  of  Reference  Benchmark 
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Procedures  for  establishing  the  fundamental  benchmark  to  serve  as  the  origin  for  the  regional 
vertical  datum  and  the  problems  interent  in  the  determination  of  local  mean  sea  level  for  that 
purpose  are  discussed  in  Vanicek  and  Krakiwsky  (1986,  Chapter  19,  p.  424-425). 

For  a  surface  point  P,  the  potential  difference  to  the  geoid  (local  vertical  reference  datum 
through  the  point  Pq  with  potential  W  =  Wq)  is, 

Cp  =  Wo  -  Wp  =  -  JdW  =  Jgdn  (2.4) 

Po  Po 

where  Cp  is  defined  as  the  geopotential  number,  which  can  be  determined  by  leveling  and 
surface  gravity  values,  starting  from  the  zero  point  Pq.  For  the  definition  of  other  variables  used 
in  equation  (2.4)  refer  to  Figure  2.3  below: 


Figure  2.3  Geoid,  Level  Surfaces,  Height  Increments  dn  and  Gravity  g 

To  compute  ^  it  is  not  necessary  to  know  the  mean  gravity  g  between  P  and  Po  along  the 
plumb  line.  Cfeopotentiai  numbers  are,  therefme,  the  only  quantities  that  can  be  determined 
without  any  density  hypothesis.  The  geopotendal  number  is  independent  of  the  particular 
leveling  line  used  for  relating  the  point  P  to  the  reference  benchmark.  The  geopotential  numter 
is  the  same  for  all  points  of  a  level  siuface  and  can  thus  be  considered  as  a  natural  measure  of 
height,  even  if  it  does  not  have  the  dimension  of  a  length  (Heiskanen  and  Moritz,  1967).  Some 
of  the  height  datums  in  the  world  today,  for  example  NAVD88,  UELN73,  etc.,  consider 
geopotendal  number  as  the  basic  quandty  for  the  dednidon  and  computadon  of  their  height 
systems.  The  geopotendal  number  C  is  measmed  in  geopotendal  units  (g.p.u.)  where 


1  g.p.u.  =  1  kgal  meter  =  1(X)0  gal  meter 


and  1  g.p.u.  coRcsponds  to  a  change  in  height  of  id)out  1  meter. 

2.1.2  OnhoBicnic  Heights 

The  orthometric  height  Hp  of  a  point  P  is  defined  as  (Heiskanen  and  MOTitz,  1967,  (4.21)) 


where  Cp  is  the  geq>otential  number  at  the  point  P  and  g  the  mean  ^yily  along  the  plumbline 
(refer  to  Figure  2.4)  between  the  surface  P  and  its  correqwnding  pout  Pq  on  the  geoid.  Since 
we  need  the  density  distribution  between  the  surface  of  the  earth  a^  the  geoid  ftn*  computing  g, 
the  practical  determination  of  this  quantity  is  inqxKsible  without  additional  assuiiq>tk)ns. 

Several  hypotheses  have  been  made  about  the  variation  of  gravity  between  the  terrain 
surface  and  the  gemd  to  approximate  this  mathematically  rigorous  quantity  g.  Depending  upon 
the  methodMgorousness  of  estimating  g,  numerous  tyjws  of  orthometric  lui^ts  are  defii^  and 
used  in  various  parts  of  the  world  today.  When  wt^mg  in  a  country  or  area  where  a  specific 
type  of  orthometric  height  is  used,  we  iwed  to  understand  its  definition  and  accuracy  limitations. 
Some  of  the  different  l^ds  of  orthometric  heights  used  in  practice  are  discussed  in  the  following 
sections. 

2.1.2.1  Rigorously  Computed  Qrthometric  Heights 

Rigorous  orthometric  heights  are  those  wha«  there  has  been  actual  attempts  to  determine 
the  actual  gravity  g,  along  the  plumb  line  from  the  surface  point  to  the  corresponding  point  on 
the  geoid,  as  closely  as  possible. 


Figure  2.4  Reference  Surfaces  Used  in  the  Reduction  of  Poincarf  and  Pney 
Frcmi  equation  (2.5)  we  know  that  the  orthometric  height  at  the  point  P  is 
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where  Cp  is  the  geopotendal  number  at  the  point  P.  Q>  »  Wo  *  Wp.  and  g  is  the  mean  gravity 
along  plumb  line  from  P  to  Fq-  The  mean  gravity  g  can  be  written  as 

g  =  ^Jg(z)<lz  (2.6) 

where  g(z)  is  the  actual  gravity  at  the  arbitrary  point  Q  which  has  the  height  z  (refer  to  Figure 
2.4). 


The  simplest  method  to  compute  an  approximation  for  the  actual  gravity  gq  at  Q  is  by  using  the 
Poincar6  and  Prey  reduction  procedure  discussed  in  Section  4.3,  Heiskanen  and  Modtz  (1967). 
This  procedure  is  based  on  the  assumption  that  the  area  around  the  gravity  station  P  is 
completely  flat  and  horizontal  (Figure  2.4)  with  the  masses  between  the  geoid  and  the  earth 
surface  having  constant  density  (2760  kgAn^).  This  flat  plate  with  thickness  equal  to  the  height 
of  the  station  is  also  called  a  Bouguer  plate.  Based  on  the  Prey  reduction  procedure  we  get  a 
value  for  the  actual  gravity  at  Q  as: 


gq  =  gp  +  0.0848  (Hp  -  Hq) 
with  g  in  gals  and  H  in  kilometers. 


(2.7) 


If  we  calculate  g  as  the  mean  of  gpivity  g^  measured  at  the  surface  point  P,  and  of  gravity  gq 
computed  at  the  corresponding  geoid  point  Q  (=  Fo,  Hq  *  0)  by  Prey  reduction  then, 

g  =  j(gp+gQ)  =  ^(gp  +  gp+0.0848Hp)  =  gp+0.0424Hp  (2.8) 


From  equations  (2.5)  and  (2.8),  we  get  the  orthometric  height  at  the  point  P  approximately  as 


Hp 


Cp 

gp+0.0424H 


(2.9) 


The  orthometric  heights,  Hp,  computed  using  ^nation  (2.9)  are  called  ‘Helmert  heights’,  with 
Hp  in  meters,  the  geopotential  number  Cp  given  in  g.p.u.,  gp  in  gals  and  the  H  is  Hp  in 
kilometers. 


Helmert  heights  computed  using  equation  (2.9)  can  be  considered  as  the  sufficient 
approximation  to  true  orthometric  heights  for  most  of  the  applications.  Sometimes  in 
mountainous  areas  and  for  highest  precision,  it  will  be  necessary  to  consider  the  variation  in 
topography  in  acWtion  to  the  infinite  Bouguer  plate.  The  heights  computed  taking  topography 
into  account,  assuming  normal  gradient  for  frM-air  reduction  and  a  constant  density  of  2670 
kg/m^  are  called  ‘Niethammer  heights’.  Detailed  discussion  on  Niethammer  heights  are  given  in 
Krakiwsky  (1%5). 

Helmert  heights  and  Niethammer  heights  are  rigorously  computed  approximations  for 
orthometric  heights  used  in  some  parts  of  the  world,  where  dense  gravi^  data  terrain  data  are 
available. 
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Another  (not  so  rigorously  computed,  but  making  use  of  geopotential  number) 
orthoaietric  height  was  proposed  by  J.  Vign^  (1954);  this  so-called  'Vignal  height’  is  used  in 
some  countries.  In  this  case,  the  mean  gravity  along  the  plumb  line  is  compute  in  two  steps 
(Krakiwsky,  1965,  pp.  96-97). 

•  The  nmnial  gravity,  is  computed  on  the  spheroid  for  the  latitude,  0p,  of  point  P. 

•  This  value  is  then  reduced  upward  to  the  mid-point  of  the  plumb  line  by  a  negative  free  air 
correction, 

Thus,  we  get  the  mean  gravity  to  be  approximately. 


5L  H 

g=Y*— f=’y.,-0.30g6-^ 

Then  the  Vignal  height  at  point  P  will  be 

u  —  Cp  _ _ Cj _ 

*  g  Y*p-0-1543  Hp 


(2.10) 


where  Y*p  is  in  gals,  Hp  is  in  km  and  Hv  in  meters. 

2. 1.2.2  Normal  Qrthometric  Heights 

Normal  orthometric  heights  can  be  considered  as  a  coarse  approximation  to  true 
Qrthometric  heights.  In  determining  this  kind  of  orthometric  height  geopotential  numbers  are  not 
used.  Most  of  the  countries  which  lack  dense  gravity  data  coverage  ^opt  this  kind  of  simple 
orthometric  height  obtained  by  the  following  relation 


Hi  =  Hi.i-t-Ah4dH  (2.11) 

where  Hi  is  the  normal  ordimnetric  height  of  the  i<^  station, 

Hi-i  is  the  nmmal  orthometric  height  of  the  previous  station  in  die  level  networit. 

Ah  is  the  leveled  difference  corrected  for  all  systematic  errors,  and 

dH  is  called  the  normal  orthotnetric  correction,  which  accounts  for  the  effect  of 

convergence  of  equipotential  surfaces  of  the  normal  field  and  is  given  by  Rapp 
(1961)  as 


dH  = -f*sin2<))Hd<|>  (2.12) 

In  equation  (2.12),  f*  is  the  gravity  flattening  defined  as  equatitm  (2.98)  in  Heiskanen  and  Moritz 
(1967).  The  negative  sign  in  the  expression  indicates  that  the  equipotential  surfaces  converge 
toward  the  poles.  Also  the  H  in  e^n.  (2.12)  refers  to  the  mean  latitude  and  elevation  between 
the  two  stations  i  tmd  (i-1),  respectively,  and  d^  is  the  difference  in  latitude  between  the  two 
stations. 
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Normal  heights  are  defined  in  the  normal  gravity  field,  that  is.  with  respect  to  a  level 
ellipsoid  with  nonnal  potential  U  =  Uo.  by  the  fcnmda 

H*  =  &,  /itfH*  (2.13) 

The  mean  nonnal  gravity  value  y  can  be  computed  without  any  hypothesis.  The  integration 
from  ‘0’  implies  that  the  normal  height  is  zero  at  the  reference  surface  from  which  the  potential 
difference,  Cp,  was  determined  by  leveling. 


Figure  2.5  The  Telluroid  (£),  Normal  ffeight  H*  and  Height  Anomaly  C 

As  depicted  in  Figure  2.5,  the  normal  height  H*  of  the  point  P  is  defined  as  the  height  of 
the  corresponding  point  Q  on  the  telluroid  (£)  above  the  reference  ellipsoid  measured  along  the 
normal  plumb  line  from  Q  toQ^o.  However,  since  the  curvature  of  the  nonnal  plumb  line  is  very 
small  and  regular  (Heiskuen  and  Mtnitz,  1967),  the  length  measured  along  the  normal  plumb 
line  can  be  considered  to  be  equal  to  the  straight  ellipsoidid  normal  hei^t  to  sufficient  accuracy. 
The  normal  height  H*  is  also  interpreted  as  the  height  of  the  surface  point  P  leading  to  the  point 
P6  on  the  so-caUed  quasi-geoid  which  is  not  a  level  surface  but  is  considered  to  be  located  close 
to  the  geoid  since  Qs  N.  Hie  telluroid  (£)  is  defined  as  the  surface,  sufficiently  close  to  the 
earth's  surface,  whose  normal  potential  U  at  every  point  Q  is  equal  to  the  actual  potential  W  at 
the  corresponc^g  point  P,  so  that  Uq  =  Wp,  corres^nding  points  P  and  Q  being  situated  on  the 
same  ellipsoidal  normal  (Heiskanen  and  Moritz,  1967).  The  height  anomaly  ^  is  the  distance  of 
the  telluroid  point  Q  form  the  earth's  surface  point  P,  or  the  height  of  the  quasi-geoid  point  P6 
from  the  reference  ellipsoid. 


The  nonnal  height  H*,  and  hence  die  telluroid  £.  can  tqipoximately  be  d^eraiined  by 
leveling  ctnnbined  with  gravity  measurements  according  to  Hetskanen  and  Mtwitz  (sec.  4-5, 
1967).  First  the  get^tential  number  of  P,  Cp  »  Wo  -  Wp  is  conqiuted  by 

P 

Cp=Jgdn  (2.14) 

Po 

where  g  is  the  locally  measured  gravity  and  dn  is  the  leveling  increment  The  ntnmal  height  H* 
is  then  related  to  Cp  by  the  analytical  expression  given  as  equation  (4-44)  (ibid.). 


H** 


Yo 


l-i-(l+f+m-2fsin*6)-^+ 

'*Yo 


(2.15) 


where  Yo  is  the  nonnal  gravity  at  the  ellipsoid,  for  the  same  ladtude  <)>.  a,  f  are  the  semi-major 
axis  and  flattening,  respectively,  for  the  ellipsoid  and  m  is  the  ctmstant  defined  by  eqn.  (2-70)  in 
Heiskanen  and  Moritz  (1967).  From  equation  (2.13)  it  is  clear  that  H*  can  be  computed 
independent  of  the  density  assumpticHi  unlike  the  o^ometric  height  H  (see  equation  (2.5)).  This 
property  of  computed  normal  heights  makes  it  preferable  to  adopt  them  as  official  heights  for 
certain  regional  vertical  networks.  For  example,  some  Eastern  European  countries,  France, 
Germany,  etc.,  use  nonnal  heights  in  their  national  height  systems. 

2. 1 .4  Inter  Comparison  Between  Various  Heights 

As  discussed  in  the  previous  two  sections,  the  basic  quantity  for  height  determination  is 
the  geopotential  number  Cp.  Depending  upon  tht  technique,  purpose  and  type  of  gravity  data 
used  we  arrive  at  different  heights.  Most  of  the  countries  in  the  world  today  adopt  orthometric 
heights  of  one  kind  or  the  other  for  their  national  height  systems.  Due  to  sparse  and  not  very 
accurate  gravity  coverage  a  few  years  ago,  nmmal  orthometric  heights  (ref.  section  2. 1.2.2)  have 
been  the  preferred  and  adopt^  heights  in  many  countries.  NGVD29  (U.S.A.),  AHD71 
(Australia),  ODN  (U.K.)  are  few  examples  of  national  height  systems  adopting  normal 
orthcKnetric  heights  in  their  regional  vertical  networks. 

Recently  a  new  vertical  datum  called  the  North  American  Vertical  Datum  of  1988 
(NAVD88)  based  on  a  minimum  constraint  adjustment  of  Canada-Mexico-US  leveling 
observations,  has  been  established  for  implementation  in  the  North  American  countries 
(Zilkowski  et  al.,  1992).  This  vertical  datum  adopts  Helmert  orthmnetric  heights  (refer  to 
Section  2. 1.2.1)  comput^  using  observed  ^vity  data  and  precise  level  differences,  which  can 
be  considered  closer  to  the  true  orthometric  heights.  Orthometric  heights  determined  using 
simple,  empirical  but  sufficiently  accurate  relations  between  geopotential  numbers  and 
orthometric  heights  are  also  used  in  countries  like  Austria  (Erker  and  Stti^el,  1989). 


Rapp  (1961)  has  carried  out  a  numerical  computation  of  orthometric  heights  at  the 
benchmarks  along  a  precise  leveling  line,  about  125  km  long,  described  as  REUN,  211 
MUnchen-Fasangarten-J-A-D-2-Mittenwald,  in  Report  no.  59  of  the  German  Geodetic 
Commission.  The  cnthometric  heights  have  been  computed  both  with  rigorous  techniques  like 
Niethammer  method  and  Helmert  method  (refer  to  Sec.  2. 1.2.1)  and  ^so  using  the  quicker 
method  of  applying  normal  orthometnc  ctmection  to  leveled  differences.  The  results  given  as 
Tables  14  and  15  of  his  report  (ibid.)  clearly  show  that  the  difference  between  the  heights 
computed  by  Niethammer  formula  and  Helmert  formula  is  on  the  order  of  1.5  cm.  But  the 
difference  between  the  Helmert  formula  and  one  obtained  by  applying  normal  orthometric 
correction  to  leveled  differences  shows  a  considerable  difference  of  about  8  cm  which  is  greater 
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than  what  could  be  allowed  on  the  basis  of  leveling  accuracy,  over  a  leveling  line  length  of  12S 
km.  This  shows  that  the  normal  orthomettic  ctHTcction  does  not  satisfactorily  take  into  account 
the  effect  of  gravity  on  the  leveling  measurements,  but  it  does  reduce  the  discrepancy  in 
elevations  due  to  non-parallelism  of  equipotential  surfaces. 

In  countries  belonging  to  former  Soviet  Union.  East  European  countries.  France. 
Germany,  etc.,  the  normal  height,  H*.  (discussed  in  Section  2.1.3)  is  used  as  official  flight  in 
their  national  height  systems.  The  relationship  between  the  normal  height  H*  and  the 
orthometxic  hei^t  H  is  ^scussed  in  Sections  (8-12)  and  (8-13)  of  Heiskanen  and  Moritz  (1967) 
and  also  in  Section  2.3  of  Rapp  and  Balasubramania  (1992). 

2.2  Establishment  of  Regional  Vertical  Datum: 

2.2J  EMStme-Proccdurc 

The  current  practice  of  establishing  a  re^onal  vertical  datum  in  a  country  or  area,  with 
some  exceptions  like  NAVD88,  UELN73,  etc.,  is  to  connect  the  precise  leveling  nets  to  several 
tide  gauges  widely  dispersed  around  the  perimeter  of  the  net  and  carry  out  an  adjusdnent  solution 
constraining  the  mean  sea  level  at  these  tide  gauges  as  zero,  with  the  implicit  assumption  that  the 
geoid  and  mean  sea  level  coincide  at  these  ti^  gauges. 

As  discussed  in  (Chapter  I  and  also  with  reference  to  Figure  2.2,  we  see  that  the  mean  sea 
levels  computed  at  the  tide  gauges  do  not  coincide  with  the  geoid  due  to  the  existence  of  sea 
surface  topography  (SST)  caused  by  ocean  currents,  temperature,  and  density  variations,  as  well 
as  air  pressure  and  wind  stress.  In  addition,  close  to  the  shore,  the  sea  bed  topography  and  river 
discharge  may  also  play  a  significant  role.  Since  the  magnitude  of  the  estimate  SST  varies  from 
one  tide  gauge  to  another,  the  resulting  adjustment  of  the  leveling  data  to  fit  these  different  mean 
sea  levels  will  lead  to  distortion  in  the  level  networks. 

To  understand  the  extent  of  variation  in  mean  sea  levels  due  to  SST  at  different  tide 
gauges,  refer  to  Figure  2.6  below  (refers  to  Figure  9,  Zilkowski,  et  al.,  1992). 


Figure  2.6  Height  Differences  Between  NAVD88  and  the  Qithometric  Height  of  Tidal 
Benchmarks  Above  LMSL  (epoch  1960-78)  (units  =  cm) 
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Figure  2.6  shows  the  height  differences  between  the  (»thotnetric  heights  of  tidal  bench 
marks  above  local  mean  sea  level  (epoch  1960-78)  at  different  tide  gauge  sites  along  the  U.S. 
Coast,  and  the  recently  made  available  NAVD88  heights.  NAVD88  vertical  datum  has  been 
established  by  a  minimum  constraint  adjustment  of  Canadian-Mexican-U.S.  leveling 
observations  holding  fixed  the  height  of  primary  benchmark  at  Father  Point/Rimouski,  Quebec. 
Canada  as  zero  height.  From  Figure  2.6,  we  can  understand  that  though  there  can  be  several 
reasons  for  the  not^  height  differences  at  the  ddal  benchmarks,  the  main  contribution  to  the 
differences  is  the  SST,  which  varies  to  some  extent  from  one  tide  gauge  to  another  along  the 
same  coast  and  vary  to  a  greater  extent  from  one  coast  to  another.  Wi^out  removing  the  SST 
effects  at  tidal  stations,  it  would  be  incorrect  to  constrain  the  heights  of  tidal  benchmarks  to  zero 
height  in  the  adjustment  solution  for  establishing  a  regional  vertical  datum. 

2.2.2  Modeling  for  Ideal  Regional  Vertical  Datum  Definition 

2.2.2. 1  Estimation  of  Sea  Surface  Topography 

One  obvious  solution  to  the  problem  of  relating  the  reference  surface  at  the  tide  gauge 
sites  to  the  mean  sea  level  is  to  model  the  SST. 

Several  global  and  regional  models  have  been  developed  so  far  mainly  using  altimetric 
data  to  estimate  the  SST.  The  Ohio  State  1991  global  SST  models  to  degree  10  and  15  have 
been  developed  based  on  analysis  of  one  year  Geosat  altimeter  data  from  the  orbits  computed 
with  GEM'T2  potential  coefficient  model  and  consistent  station  coordinates  (Rapp  et  al.,  1991). 
The  commission  error  in  the  estimation  of  SST  in  the  ‘ocean  areas’  is  about  ±10  cm  in  the  case 
of  10,10  model  while  the  error  for  the  15,15  solution  is  higher  because  twice  as  many 
coefficients  are  included  in  the  degree  15  model.  The  extrapolation  of  SST  at  the  tide  gauge 
points  from  these  models  is  not  recommended  at  this  time  since  the  inaccuracy  in  the  estimate  of 
SST  in  shallow  coastal  areas  will  be  too  high  to  contribute  anything  to  adjusted  level  heights. 

Efforts  have  been  made  by  icw  researchers  to  model  the  local  SST  variation  by 
considering  monthly  mean  sea  level  and  other  data  for  several  tide  gauges.  Merry  and  Vanicek 
(1983)  proposed  a  technique  which  makes  use  of  the  response  of  sea  level  to  external  effects 
such  as  atmospheric  pressure,  temperature,  wind,  and  so  on.  Their  experiment  with  relevant  data 
at  various  tide  gauge  sites  in  the  Maritime  Provinces  of  Canada  concluded  that  their  results  do 
not  represent  the  entire  SST  at  the  tide  gauges  but  only  reflect  the  local  vturiations  of  SST. 
Similar  studies  have  been  carried  out  by  Merry  (1990)  to  determine  the  variation  in  mean  sea 
levels  at  the  tide  gauge  sites  along  the  coast  of  South  Africa.  Heck  and  Rummel  (1990)  have 
developed  formulas  for  extrapolating  current  and  SST  information  along  shallow  coastal  waters 
from  locally  measured  data.  The  main  difficulties  with  these  methods  are: 

•  They  are  better  suited  to  the  determination  of  SST  differences  between  adjacent  tide  gauges 
rather  than  to  the  determination  of  absolute  SST. 

•  Due  to  the  inadequacy  of  the  hypothesis  and  imperfect  modeling  of  frictional  forces  in  the 
continental  shelf  areas  (eg.  bottom  friction  and  surface  windstress)  the  SST  required  for 
relating  vertical  reference  surface  to  the  tide  gauge  sites  cannot  be  estimated  to  intended 
precision  (Heck  and  Rummel,  1990). 

•  The  effort  required  to  obtain  the  necessary  data  is  significant  and  the  methods  may  prove  to 
be  non-viable  purely  on  economic  grounds. 
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2.2.2.2  Combined  Adjustment  Solution  for  Defining  Regional  Vertical  Datum 

Several  techniques  for  defining  a  regional  vertical  datum  have  been  put  forward  by 
various  authors  in  recent  years.  The  general  consensus  amongst  them  is  that,  for  establishing  a 
regional  vertical  datum,  we  need  to  consider  the  potential  differences  derived  f^m  a  combination 
of  leveling  and  gravity  data,  precise  altimeter  data  from  missions  like  ERS-l, 
TOPEil^POSEn[)ON,  etc.,  and  the  tide  gauge  records  at  the  primary  tide  gauge  stations  and  come 
up  with  a  combined  adjustment  solution. 

Groten  and  Muller  (1990)  have  suggested  such  a  combined  adjustment  model  using 
gravity  field  as  well  as  leveling  data,  altimetry  and  tide  gauge  data  and  three-  dimensional 
coordinates  directly  determined  by  GPS.  The  ongoing  Bass  Strait  vertical  datum  definition 
project  in  Australia  (Rizos  et  al.,  1991)  is  also  another  effort  to  integrate  geometric  information 
(from  space  geodetic  techniques)  and  gravity  field  information  (from  surface  ^vity  data  and 
accurate  regional  geoid  models)  for  regional  vertical  datum  definition.  Vanicek  (1991)  has 
proposed  five  different  approaches  for  defining  a  regional  venical  datum  depending  on 
availability  of  data  and  rigorousness  required  in  defining  the  vertical  datum. 

Based  on  the  assumption  that  realistic  models  to  estimate  the  SST  at  reference  tide 
gauges  may  be  available  in  ^e  near  future,  the  following  procedure  can  be  adopted  to  come  up 
with  an  adjustment  solution  for  defining  the  vertical  datum.  Assuming  that  we  have  the 
following  data: 

•  observed  geopotential  differences  between  the  leveling  bench  marks  expressed  in 
geopotential  units  (gpu),  ACy  =  Cj  -  Ci  for  all  i,  j  (except  i  =  j) 

•  geopotential  numbers  at  the  reference  tide  gauge  sites  obtained  from  realistic  models  to 
estimate  the  sea  surface  topography. 

We  can  set  up  an  observation  ^nation  model,  introducing  the  estimated  geopotential  numbers  at 
the  tide  gauge  sites  as  stochastic  constraints  The  formulation  of  the  model  will  be 


y=AX-»-e,D{y}  =  a5P-‘ 

Zq  =  KX+eo,D|zoJ  =  OoPo* 


(2.17) 


where  y 
A 

m 

X 


Zo 


Co 

p 

Po 

K 


D 


is  the  nxl  observation  vector  (of  geopotential  differences  ACy) 

is  the  nxm  coefficient  matrix  of  observations 

rank  A  =  q^m^n 

is  the  number  of  parameters 

is  the  mxl  parameter  vector  (unknown) 

is  the  nxl  random  error  vector 

is  the  £xl  vector  of  stochastic  constraints  having  the  same  structure  of  observation 
vector  with 

as  corresponding  ^xl  random  error  vector 

is  the  weight  matrix  of  observations 

is  the  weight  matrix  of  the  stochastic  constraints 

is  the  common  variance  component  (unknown) 
is  the  txm  coefficient  matrix  of  the  constraints 


rank  K  =  /  <  m  and  rank 


A 

K 


=  iti 


is  the  dispersion  operator. 
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Assuming  no  correlation  between  y  and  Zq.  we  can  write  the  least  squares  solution  for  the 
adjusted  parameters  to  be  (B.  Schafhm,  1990) 


X  =  (A^PA+KTPoK)‘’(ATpy+KTPoZo)  (2.18) 

and  the  dispersion  matrix  of  the  estimated  parameters  to  be 

D{x}=a2(A^PA+K’'PoK)''  (2.19) 

The  problems  with  solution  (2.18)  at  present  conditions  as  pointed  out  by  Zilkoski  et  al. 
(1992)  are: 

•  Realistic  models  to  compute  the  SST  at  the  tide  gauge  sites  to  required  accuracy  are 
unavailable  at  this  time,  and 

•  Ana  priori  estimate  of  standard  errors  for  tidal  stations  where  SST  effects  were  not  removed 
completely  would  be  too  large  relative  to  precise  leveling  differences  to  allow  any 
contribution  to  the  final  adjusted  heights. 

To  overcome  the  inconsistencies  in  the  estimated  SST  values,  the  ‘robust  sequential 
updating  procedure’  suggested  by  B.  Schaffrin  (1987,  1989,  1990)  can  provide  an  excellent 
adjustment  technique  in  the  presence  of  ‘weak  prior  information*.  Hiis  procedure  provides 
solutions  when  the  stochastic  prior  information,  included  as  stochastic  constraints  in  our 
observation  equation  model,  is  allowed  to  have  a  certain  (guessed)  bias  Po  or  even  certain 
suspected  scaling  error  in  the  bias  po.  As  discussed  in  B.  Schaffrin  (1987),  the  adjusted  solution 
suppresses  the  inconsistencies  in  the  estimated  geopotential  numbers  at  the  tide  gauges  without 
destroying  the  homogeneity  of  the  geopotential  differences  observed  at  the  leveling  benchmarks. 

Vanicek  (1991)  has  suggested  using  orthometric  heights  computed  with  equation  (2.1)  at 
space  geodetic  stations  can  be  included  with  leveling  data  in  a  leveling  network  adjustment,  in 
o^er  to  homogenize  the  error  distribution  within  the  network.  But  since  the  uncertainties  of 
geoid  height  differences  used  to  convert  the  ellipsoidal  height  differences  to  orthometric  height 
differences  are  presently  large  compared  to  fonn^  errors  of  leveling  height  differences  (Zilkoski 
et  al.,  1992),  such  a  suggestion  is  impracticable  as  of  today. 

2.2.2.3  Summary  and  Ctonclusions 

1.  The  primary  reference  surface  for  heights  is  the  geoid.  Since  it  is  approximated  as  the  mean 
sea  level  (that  is,  as  a  temporal  mean  of  all  the  instantaneous  sea  level  heights  at  the  epoch 
‘t’),  it  is  imperative  that  the  tide  gauge  records  be  included  in  the  regional  vertical  danim 
definition. 

2.  The  mean  sea  level  departs  from  the  geoid  due  to  the  existence  of  the  quantity  called  Sea 
Surface  Topography  (SST).  Since  this  quantity  dejwnds  on  meteorological,  geophysical, 
hydrological  and  oceanographic  components  at  the  tide  gauge  sites  and  its  magnitude  is  a 
combination  of  all  these  components,  it  varies  considerably  from  one  region  to  another. 
Therefore,  it  is  impracticable  to  include  SST  for  a  Global  vertical  datum  modeling.  But  it  is 
possible  that  the  change  in  SST  from  one  tide  gauge  location  to  anoAer  in  a  given  datum  can 
be  determined  to  a  reasonable  accuracy,  which  could  enable  multiple  relative  tide  gauge 
information  to  be  incoiporated  into  a  regional  venical  datum  definition  without  creating  the 
distortion  problem  noted  earlier. 
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3.  At  present  suitable  techniques/realistic  models  to  determine  the  SST  exactly  to  intended 
precision  are  not  available.  But  without  removing  SST  el^ects  at  the  tidal  stations,  it  would 
also  be  incorrect  to  constrain  the  heights  of  tidal  benchmarks  to  zero  height  in  the 
establishment  of  regional  vertical  datums.  In  such  a  scenario,  the  suggested  approach  will  be 
to  estimate  the  SST  as  accurately  as  possible  and  then  set  the  tidal  tenchim^  height  to  be 
equal  to  (MSL-SST).  Now  the  tide  gauge  heights  (or  geopottential  number)  can  be  used  as 
stochastic  constraints,  keeping  in  mind  the  noise  at  these  stations  will  be  too  high  due  to  the 
undetermined  bias  in  the  SST  values.  To  accommodate  for  the  bias  in  the  estimated  SST 
values  at  the  tide  gauge  sites  and  possible  scaling  error  also  in  the  bias,  the  robust  sequential 
updating  procedure  suggested  by  B.  Schaffrin  (1987)  can  be  adopted  for  adjustment  of 
regional  vertical  datum. 

4.  The  suggestion  to  include  the  orthometric  heights,  computed  using  space  geodetic  techniques 
and  global/regional  geopotential  models,  in  order  to  homogenize  the  error  distribution  within 
the  regional  vertical  datum  adjustment  seem  to  be  impracticable  as  of  today,  as  the 
uncertainties  of  geoid  heights  (or  geoid  height  differences)  are  presently  large  compared  to 
formal  errors  of  leveling  height  differences  (Zilkoski  et  al.,  1992). 
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CHAPTERS 


GLOBAL  VERTICAL  DATUM  DEFINITION « AN  IDEAL  APPROACH 


Defining  and  realizing  a  Global  Vertical  Datum  to  an  accuracy  of  ±10  cm  is  essential 
today  for  the  various  reasons  cited  in  Chiqiter  L  With  the  existing  measuring  techniques  and 
available  improved  models  and  procedures,  the  question  is  not  whether  we  n^  to  establish  a 
precise  vertical  datum  that  can  be  used  <mi  a  global  basis,  but  to  what  accuracy  uul  how  soon. 

Recognizing  the  fact  that  the  data  are  neitho*  uniformly  available  all  over  the  world  nw  to 
the  required  precision,  noodelling  procedures  fcv  a  Global  Vertical  Datum  that  can  be  realized  are 
discussed  in  this  chiqiter.  The  available  data  as  of  today  can  be  used  in  the  models  developed  for 
estimating  the  parameters  that  will  define  the  Global  Vertical  Datum,  but  only  as  a  tirst  ito^tion 
attempt.  As  understandable,  it  may  require  several  iterations  before  achieving  an  Ideal  Global 
Vertical  Datum  that  could  meet  all  die  requirements  stated  in  Chapter  1 . 

The  modelling  procedures  have  been  developed  for  different  data  scenarios  but  with  the 
assumption  that  the  data  are  available  everywhere  and  also  to  the  required  precision,  that  a 
Global  Vertical  Datum  established  to  an  accuracy  of  ±10  cm  would  warrant  It  is  also  assumed 
that  the  regional  vertical  datums  have  been  established  based  on  the  procedures  discussed  in 
Chapter  2. 

1  Data  Reouiiements  for  ModeUing  Purposes 

Basically  four  kinds  of  data  are  essential  for  noodelling  the  Global  Vertical  Datum.  They 
are  accurately  determined  free  air  gravity  anonudies,  precise  heights  of  the  stations  above  the 
regional  vertical  reference  datums,  an  accurate  global  geopotential  model,  and  accurate 
elUpsoidal  heights  of  the  stations  above  an  adopted  geocentric  ellipsoid.  In  this  section  a  full 
description  and  accuracy  requirements  of  requirra  data  types  are  given. 

3.1.1  Free-air  Gravity  Anomalies 

Two  different  approaches  are  used  in  defining  the  free-air  gravity  nomalies.  The  first 
approach  is  termed  as  'a  classical  qipoach'  (Chapter  2-13,  Heiskanen  and  Moritz,  1967)  and  the 
other  as  Moloctensky  approach'  (Section  42,  Moritz,  1980). 
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Figure  3.1  Various  Equipotendal  Surfaces  and  Heights  Used  in  Defining 
the  Free-air  Gravity  Anomalies 

In  the  classical  approach  the  gravity  anomaly  Age  is  defined  as  the  difference  in 
magnitude  of  the  gravity  vector  at  Pq  and  the  normal  gravity  vector  at  (refer  to  Figure  3.1). 

(3.1) 

where  g«  is  the  magnitude  of  the  observed  gravity  at  P  reduced  to  geoid  point  Pq,  and  y.  is  the 
magnituw  of  noirntd  gravity  at  the  point  Qq  cmnputed  on  the  reference  ellipsoid. 

The  observed  gravity  at  the  point  P  is  reduced  to  the  point  Pq  on  the  geoid  using  free-air 
reduction  procedures  discussed  in  Heiskanen  and  Moritz  (1967,  Section  2.23).  If  gp  is  the 
observed  ^vity  at  the  surface  of  the  earth  then  the  value  of  gp  at  the  geoid  may  be  obtained  as 
a  Taylor  expansion:  ** 


8g 


H-. 


(3.2) 


ds 

where  ^  is  the  vertical  gradient  of  the  gravity  at  P  and  H  is  the  orthometric  height  of  the  station 

above  the  geoid.  Considering  only  the  linear  term,  the  quantity  accounts  for  the  free- 

air  reduction  of  the  observed  gravity  to  the  geoid.  Here  we  have  assumed  that  there  are  no 
masses  above  the  geoid  or  that  such  masses  have  been  removed  beforehand  so  that  this  reduction 
is  indeed  carried  out  in  "free-air".  From  equations  (2-215)  and  (2.217)  in  Heiskanen  and  Moritz 
(1967),  we  note  that  the  vertical  gradient  of  gravity  can  be  computed  by: 


dg  _  dy  ^  3Ag 

5h  ah  ah 


(3.3) 
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where  ^  is  the  normal  gravity  gradient  given  by  equation  (2*121)  and  is  the  anomalous 


part  which  can  be  computed  using  equation  (2-217)  in  Heiskanen  and  Moritz  (1967).  Due  to 
lack  of  dense  gravity  coverage  required  for  cmnputing  the  actual  vertical  gra^ent  of  gravity 
using  equation  (3.3),  it  is  generally  approximate  by  normal  gradient  of  gravity,  setting  the 
anomalous  part  equal  to  zero.  Ba^d  on  this  apim>ximation,  equation  (3.2)  can  written  to 
second-order  as: 


Substituting  the  values  for  the  first  and  second  derivative  of  normal  gravity  from  equations  (2- 
121)  and  (2-122)  from  Heiskanen  and  Moritz  (1967),  in  equation  (3.4),  we  can  write  the  free-air 
gravity  anomaly  in  'classical  sense*  as; 


Ag.  =  gp  -  f  1  -  2{H-  f  +  m  -  2f  sin^  «|))M + j' 


(3.5) 


In  the  Molodensky  approach,  the  gravity  anomaly  is  defined  as  (refer  to  Figure  3.1): 


AgM  =  8p  ■  Yq  (3.6) 

where  gp  is  the  ^vity  observed  at  the  surface  point  P  and  the  Yq  is  the  normal  gravity  at  the 
corresponding  point  Q  on  the  telluroid  (2).  The  normal  gravity  at  the  point  Q,  Yq,  is  computed 
from  the  normal  gravity  at  the  ellipsoid  by  the  normal  free-air  reduction,  hm  now  applied 
upwards. 


dy  *  1  9*y 

YQ=YQo+-g^H 


(3.7) 


where  H*  is  the  normal  height  of  the  point  P.  A  direct  formula  for  computing  Yq  at  Q  is  also 
given  by  equation  (2-123)  in  Heiskanen  and  Moritz  (1967)  as; 


Yq  =  Yqo 
From  (3.6)  and  (3.8)  we  write: 


1  -  2{H- f -I- m  -  2f  sin*  3|^^  J 


(3.8) 


AgM=gp-YQ„  l-2(l+f+m-2fsin*<|»)^+3^^J 


(3.9) 


A  reasonable  estimate  for  the  numerical  difference  between  the  free-air  anomalies  computed 
using  equations  (3.5)  and  (3.9)  can  be  determined  as  follows; 


Considering  only  up  to  the  first  order  normal  gravity  gradient, 
can  write  from  equations  (3.7)  and  (3.4)  as: 


dh 


=  0.3086mgal  /  meter ,  we 
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(3.10) 


“Yob  '♦■0.3086H 
Age  =gp  -Yqb  +0.3086H 

AgM-Ag.  =0.3086(H*-H)  (withyQ,  ^Yg,) 

= 0.3086(h -C  -  h  +  N)  =  0.3086(N  -  Omgals 

where  N  is  the  geoid  undulation  and  ^  the  height  anomaly.  Prom  equation  (8-104)  in  Heiskanen 
and  Moritz  (1967),  we  get: 

(C-N),««=-K).lH£.Hta  (3.11) 

where  H  is  the  elevation  of  the  station  above  vodcal  reference  datum  and  H^y  the  average  height 
of  the  area  ctmsidered. 

For  an  elevation  of  1000  meters,  substituting  (3.1 1)  in  (3.10)  we  get: 

AgM  *  Age  ~  -0.03086  mgals  (3. 12) 

Though  the  numerical  difference  between  the  two  anomalies  as  implied  by  equation 
(3.12)  is  sm^l,  conceptually  they  are  very  different.  While  the  free-air  anomalies  computed 
using  equatitm  (3.9)  are  referred  to  ground  level,  die  classical  free-air  anomalies  from  (3.5)  are 
referred  to  sea  level  (Heiskanen  and  Moritz,  1967,  p.  293).  This  distinction  should  be  carefully 
kept  in  mind  when  using  the  different  types  of  ftM-air  anomalies  in  the  definition  of  globtd 
vertical  datum. 

The  free-air  anomalies  detined  and  computed  using  either  approach  need  to  be  corrected 
for  the  following  systematic  corrections: 

(a)  Atmospheric  Conwtion  (figAli 

In  the  problem  of  physical  geodesy,  whether  following  tire  Molodensky  approach  or 
classical  approvh,  it  is  assum^  that  the  anomalous  potential  T  outside  the  boundi^  surface  is  a 
harmonic  function  or,  in  other  terms,  that  the  space  outside  the  surface  is  empty.  Therefore  the 
effect  of  atmosphere  on  Ag  must  be  removed  by  cmnputation.  Detailed  ctiscussicms  on  the  effect 
of  atmosphere  on  gravity  anomalies  and  the  amount  of  correction  needed  to  be  applied  on 
observed  gravity  can  be  found  in  Moritz  (1980,  pp.  422-425),  Pavlis  (1988, 1991)  and  Rapp  and 
Pavlis  (1990). 

We  can  make  use  of  the  following  relation  given  by  Pavlis  (1991,  p.  33)  for  computing 
the  atmospheric  cmrection  that  needs  to  be  added  to  all  gravity  anomalies  derived  tirom  terrestrid 
observations. 

SgA  =0.8658  -9.727xl0-*Hp  +  3.482xl0-’Hj(mgal)  (3.13) 

where  Hp  is  the  orthometric  height  of  the  gravity  station  in  meters. 

Since  the  indirect  effect  of  atmospheric  cOTtections  on  the  computed  height  anomaly/geoid 
undulation  is  in  the  order  of  -0.7  cm  (Moritz,  1980,  p.  425),  it  can  be  safely  neglected. 
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In  the  computation  of  gravity  anomalies,  the  normal  gravity  calculations  require  the 
values  of  geodetic  Latitude  for  tte  points  of  gravity  observation.  Since  in  the  gravity  formula  the 
latitude  values  used  are  defined  in  the  local  horizontal  datum  (which  is  not  generally  consistent 
with  global  geocentric  datum),  a  systematic  error  is  introduced  in  the  norm^  gravity  computed 
and  thus  in  the  gravity  anomaly.  As  the  transformation  parameters  from  regitmal  horizontal 
datum  to  global  geocentric  datum  can  be  determined  to  M.l  arc  second  accuracy  in  geodetic 
latitude  and  longitude  (refer  to  Tables  9.1  and  9.2,  DMA  Tech.  Report,  Part  I,  Dm.  1987),  the 
error  in  gravity  anomalies  due  to  horizontal  datum  inconsistancy  can  be  easily  computed  and 
applied.  A  detailed  discussion  on  the  effect  of  horizontal  datum  inconsistency  on  observed 
g^vity  anomalies  and  hence  on  undulations  computed  using  them  are  discussed  in  Chapter  4. 
The  magnitude  of  the  correction  to  gravity  anomalies,  5gH.  is  about  0.3  mgals  (Heck,  1990, 
Table  4)  and  the  correction  to  the  corresponding  undulation  computed  using  them  is  about  10 
cm. 

(c)  Gravity  Formula  Correction  (Scq): 

The  computed  free-air  anomalies  in  different  parts  of  the  world  should  correspond  to  the 
same  normal  gravity  formula  used  for  computing  the  normal  gravity  values.  If  they  are 
computed  with  different  gravity  formulas  then  corrections  are  needed  to  be  applied  to  make  them 
all  compatible  to  a  single  reference  formula.  For  example,  if  some  of  them  referred  to  GRS67 
formula  and  the  rest  to  GRS80  formula,  then  we  can  use  the  following  relation  given  by  Moritz 
(1992); 

(Y1980  -  Y1967)  =  (0.8316  +0.782  sin2<|)  -  0.(XX)7  sin^)  mgal  (3. 14) 

to  compute  the  correction  to  be  applied  to  the  gravior  anomalies  computed  with  GRS67  formula. 
More  discussions  on  the  correction  to  be  applM  to  currently  available  fiee-air  anomaly  data  sets 
are  included  in  Chapter  4. 

(d)  EUipsaidal 

In  most  of  the  formulas  in  physical  geodesy  relating  to  gravity  anomalies  and  geoidal 
heights  such  as  Stokes'  formula,  Molodensky  series  etc.,  one  usually  neglects  the  flattening  of  the 
reference  ellipsoid,  arriving  at  expressions  which  are  framally  valid  on  the  sphere  only.  Since 
this  spherical  approximation  can  cause  an  error  of  about  ±20  cm  in  the  computed  undulation  and 
in  our  problem  of  global  vertical  datum  definition  we  need  higher  accuracies,  ellipsoidal 
corrections  to  free-air  gravity  anomalies  should  be  applied. 

Rapp  and  Pavlis  (1990)  have  given  the  necessary  analytical  expressions  for  the 
computation  of  ellipsoidal  correction  terms  that  are  to  be  applied  to  the  free-air  gravity 
anomalies.  The  above  terms  are  of  the  order  of  tens  of  microgais  and  of  long  wavelength  in 
nature.  Their  effects  have  been  studied  in  detail  by  Cruz  (1986)  ^  Pavlis  (1988, 1991). 

1.2  Precise  Heights  Above  Regional  Vertical  Reference  Datum 


The  difierent  types  of  heights  that  are  used  in  the  various  regional  vertical  networks  in  the 
world  are  discussed  in  C^hapter  2,  the  most  common  heights  being  the  normal  heights,  Helmen 
orthometric  heights  and  normal  orthometric  heights.  As  discussed  in  Chapter  2,  the  basic 
quantity  for  height  determination  is  the  geopoiential  number  C  which  is  given  by  the  relation; 


(3.15) 


C  =  IgjAni 


where  gi  is  the  measured  gravity  and  Ani  the  leveling  increments.  In  this  nxxlelling  effort,  it  is 
assumed  that  the  leveling  increments  meet  the  first  order  leveling  procedure  requirements  and  the 
gravity  observaticms  refer  to  the  International  Gravity  Standardization  Netwo^  1971  (IGSN71) 
which  defines  die  primary  global  gravity  datum. 

3.1.3  Ellipsoidal  Heights 

Precise  ellipsoidal  heights  at  the  fundamental  stations  used  in  the  modelling  are  obtained 
from  observations  with  currently  available  space  geodetic  techniques  like  SLR,  VLB  I,  GPS  and 
DORIS  netwmics.  The  original  results  from  these  observations  will  be  the  cartesian  coordinates 
X,  Y,  Z  of  the  stations  with  respect  to  a  specific  coordinate  system.  Considering  a  global 
geocentric  reference  system  such  as  International  Terrestrial  Reference  Frame  1990  (TTRF  1990) 
as  a  basic  reference  system,  the  individual  station  coordinates  in  their  respective  coordinate 
frame  can  be  transfonned  to  the  global  reference  frame.  Boucher  and  Altamimi  (1992)  give  the 
seven  parameters  necessary  to  convert  between  some  general  reference  systems  and  ITRF90. 
For  more  discussions  refer  to  Rapp  and  Balasubramania  (1992). 

3.1.4  Global  Geopotential  Models 

The  gravitational  potential  Held  at  any  point  outside  the  earth  can  be  expressed  as  an 
infinite  series  of  spherical  harmonics  as  follows: 

VE(r,e,A,)  =  ^fl+ i  i  (c„„cosinX+s„„sinmX)P„(cos0)l  (3.16) 

*  ii*2m-0Vry  J 

where  GM  =  product  of  the  gravitational  constant  and  mass  of  the  earth, 
a  =  mean  equatorial  radius  of  the  earth, 

=  fully  normalized  spherical  harmonic  coefficients. 

Pan  (cos6)  =  fully  normalized  Legendre  polynomial, 

n,  m  =  degree  and  order  of  spherical  hamwnic  expansion,  and 

r,  9,  A.  =  spherical  coordinates. 

In  practice  since  the  gravitational  potential  has  to  be  determined  from  the  combination  of  satellite 
and  terrestrial  information,  we  can  represent  the  potential  only  by  a  truncated  series  with  a  finite 
number  of  terms.  Such  a  Hnite  expansion  of  potential  coefficients,  complete  to  degree  and  order 
Nmax,  will  be  called  as  a  global  geopotential  model  and  equation  (3.16)  becomes: 

V(r,e,X)=^^l+”r  iJiJ(c„cosmX+LsinmX)p„„(cose)j  (3.17) 

The  higher  the  maximum  degree,  Nmax>  in  the  series,  the  greater  the  resolution  we  have  in 
representing  the  potential  and  quantities  derived  from  it.  For  a  series  that  is  complete  (i.e.,  all 
coefficients  are  Imown)  to  degree  N,  the  resolution  (R)  of  the  geoid  model  can  be  estimated  from 
(Rapp,  1991): 


R(km)  = 


lllkmxl80_19980(km) 
N  N 


(3.18) 
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Pertinent  information  cm  available  ^obal  ge(^x}tential  models  are  given  in  Featherstone  (1992, 
Table  3.1),  and  they  ate  reproduced  in  Table  3.1  below. 

Table  3.1:  Information  on  Global  Geopotential  Models. 


ftoW  ftcpoltniial  aodtt 

B9 

CU^  (.»|**) 

(■) 

data  aaad 

0SU91A  (lUpp  tt  of.  1991) 

360 

3986004J6X10* 

6378137 

«i*».»8 

0SUS9A  (Kapp  It  PavUt,  1990) 

360 

3986004 J8x  10* 

6378137 

OSU89B  (lUpp  It  PavlU,  Aid) 

360 

3986004  J6x  10* 

6378137 

«4A*.tS48 

OSUaSE  (Rapp  It  Cnx.  1986) 

360 

3986004J6X10* 

6378138 

OSU86P  (Rapp  k  Cnu.  IW 

360 

3986004 J6x  10* 

6378138 

«4AX,tS4S 

OSUSl  (Rapp,  1981c) 

160 

3986004  J6x  10* 

6378137 

tt4a,lc 

GEM-Tl  (Manh  <t  ol.  1987) 

36 

3986004  J6X 10* 

6378137 

St 

GEM-T3  (Manb  <i  «/.  1990a/b) 

36 

398600446x10* 

6378137 

tt 

GEM-T3  (Nortm  «  «/.  1991) 

30 

398600446x10* 

6378137 

GRIM3  L1  (Roisbor  c<  «(.  1983) 

36 

3986003x10* 

6378140 

*tA».tS 

CRl^U-Sl  (Rctcbor  k  Balmiao,  1991) 

36 

3966004.4x10* 

6378136 

WGS  St-EGM  (Komar,  1984) 

IS 

3986003x10* 

6378137 

«t4a,tf 

GPM2  t  (WcQXci.  1983) 

300 

'  — 

— 

IFESSEl  t<  (BtaV  «<  </.  1989) 

360 

— 

— 

*tAX.tS 

GRS  80  (Moriu,  1980a) 

— 

3986003x10* 

6378137 

— 

ii  uU«  3.1:  St  -*  sMlIiU  orbit  poriurbuioa  iatemiKm, 

M  —  uttOit*  iltiaottr  4tfivt4  iafematioB, 
t(  —  obstr««4  tartMtiUi  grwitx  aLBOouiits.  wd 
pf  —  pitdicud  tonoftriai  (rovitjr  aaomSoo. 
t  aot  cwTtatlr  bold,  t  tailoiod  oMdtL 

3.2  Modelling  for  Vertical  Datum  Definition 

3.2.1  Procedure  when  Free-air  Anomalies  in  Motodenskv  Sense  and  Normal  Heights  of 
Stations  are  Known 

The  following  procedure  can  be  adopted  in  setting  up  the  observation  equations  for 
defining  the  Global  Vertical  Datum  when  normal  heights  of  the  stations  referred  to  regional 
vertical  datum  and  the  free-air  anomalies  in  the  Molodensky  sense  are  known; 


Global  Vertical  Datum 

w  =  wA/ 


•w=w; 


Regional  Vertical  Datum 


Figure  3.2  Various  Reference  Surfaces  and  Normal  Height 

In  the  above  figwe,  Qi  is  the  fundamental  benchmark  serving  as  the  ori^n  for  the 
regional  vertical  datum  'i*.  The  prcKedure  for  establishment  of  this  benchmark  is  discussed  in 
Chapter  2.  Let  be  the  nonnal  height  of  the  point  P  coirqiuted  with  geopotential  number  Cp  = 
Wo  -  Wp,  where  Wq  is  the  geopotf'ntial  of  the  unique  (but  unknown)  fundamental  surface 
referred  to  as  the  Global  Vertical  Datum,  and  H*  is  the  nonnal  height  of  the  same  point  P,  but 
computed  with  geopotential  number  Cp*W‘o-Wp,  referring  to  ith  regional  vertical  datum. 
From  equation  (2.15): 


(3.19) 


r*  c'  f  c' 

h;=-^  l+(l+f+m-2fsin*<^)-^+ 

Yqo  ^Qo 

i+(i+f+in-2fsin*<(>)-^+f-^1 
Yqo  '  (aYo^J 


andH:-H;-f^-^  =^(Cp-Cp) 

where  Yqo  is  the  normal  gravity  at  the  ellipsoid,  for  the  same  latitude 
Denoting  Cp  -  Cp  =  W,,  -  Wj  =  Cq.o  ,  then 


y 


(3.20) 
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From  equation  (8-8)  in  Heiskanen  and  Moritz  (1967),  the  normal  ^vity  for  die  lelluroid 
point  Q  is  given  by 


YQ=Yg,+ 

Neglecting  higher  order  terms, 

ei-iOL'*"'-' 


-Yq  + 


(3.21) 


r»yl  .fa'rl 

fHf  +  HM 

Uh-J 

0.1  2  J 

(3.22) 


where  y^Ms  the  normal  gravity  computed  ftom  elevation  information  referring  to  the  regional 
vertical  oatum  (i.e.  Hj). 

From  the  previous  section  we  know  the  fiee-air  anomaly  in  Molodensky  sense,  Agj^,  is 
defined  as: 


AgM  =  8p  *  Yq  (fro®  equation  (3.6)) 

where  gp  is  the  magnitude  of  the  actual  gravity  acceleration  at  the  surface  point  P,  and 

Yq  is  the  magnitude  of  the  normal  gravity  aa:eleration  at  the  corresponding  telluroid  point 

Q. 

Therefore,  we  can  write  from  equations  (3.20)  and  (3.22) 


Ag»,=Agg,)-^ 

Yqo 


(3.23) 


where  AgM.Agjjf  are  the  gravity  anomalies  cmnputed  with  normal  gravity  yq  determined 
using  normal  heists  Hp.H*  of  the  point  P  refierring  to  Global  Vertical  Datum  and 
regional  vertical  datum  ’i’,  respectively. 

Cq,o  is  the  potential  difference  between  the  geopotential  surfaces  referred  to  as  global 
and  regional  vertical  datums. 


Now  define 


(3.24) 
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From  equations  (2-121)  and  (2-122)  in  Heiskanen  and  Moritz  (1967)  with  spherical 

2 

=  and. 


approximation  we  have, 


6 

w 


h;>H: 


(it  is  small  and  can  be  neglected). 


Therefore,  we  can  write  equation  (3.24)  to  be. 


2 

**  R 


(3.24’) 


From  equation  (43.42)  in  Moritz  (1980),  the  computation  of  the  height  anomaly  C  can  be  done 
using  the  following  expression: 

Ag„S(v)da-»-  z  g.S(Hr)do  (3.25) 

where  Co  =  ^  (refen  to  equation  (2-181)  in  ffeiskanen  and  Moritz  (1967)), 

Ry  y 

S(GM)  =  GM  -  (GM)o  is  the  difference  between  the  true  geocentric  gravitational  constant 
and  that  of  adt^ted  reference  ellipsoid  and  AWq  =  Wq  -  Uo, 

R  is  the  mean  radius  of  earth, 

a  is  the  spherical  surface, 

y  is  the  mean  gravity  of  the  earth, 

Agj^  are  the  free-air  gravity  anomalies  defined  in  Molodensky  sense, 
y  is  the  spherical  distance  between  the  computation  point  P  and  the  variable  points  K, 

S  (\|f)  is  the  Stokes'  function, 

gn  are  the  correction  terms  evaluated  by  the  recursive  formula. 

g,=-iz%(g,_.)  (3.26) 

r-J 

In  equation  (3.26),  z  is  the  hei^t  difference  between  the  station  point  P  and  varying 
points  K,  z  =  hK  -  hp.  The  operator  Ln  is  evaluated  recursively  as  (refers  to  equation  (45-47)  in 
Moritz  (1980)) 

L„(Ag)  =  lL,[L„.,(Ag)]  (3.27) 

In  planar  approximation,  Li  is  given  by 
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(3^8) 


with  /o  =  2Rsin^ 

Starting  from  gg  »  Ag|^  from  equation  (3.26)  we  can  write: 

gj  =  -  (hj^.  -  hp)  L|^(g^ 

gj  =  -  (h^  -  hp)  Li(g,)  -  (h^  -  hp)2  L2(gQ)  (3.29) 


and  so  on.  The  terms  gi,  g2.  etc.  are  called  as  Molodensky  correction  terms. 

The  evaluation  of  the  equations  (3.26)  to  (3.29)  by  formulating  convolution  integrals  for  the 
planar  approximation  and  then  solving  them  in  frequency  domain  using  Fast  Fourier  Transform 
techniques  are  discussed  in  Sideris  and  Schwaix  (1986, 1988). 

The  free-air  gravity  anomalies  derived  in  Molodensky  sense  and  used  in  equation  (3.2S)  are  to  be 
corrected  for  atmospheric  and  other  corrections  listed  in  Section  3.1.1.  We  can  denote  these 
corrected  anomalies  as  AgM. 

Based  on  the  discussions  in  Sideris  and  Schwarz  (1988),  for  a  highly  accurate 
computation  of  the  height  anomaly  using  equation  (3.2S),  consideration  of  Molodensky's 
correction  terms  up  to  g2  seems  to  be  significant.  But  numerical  tests  (ibid.)  also  show  that  these 
correction  terms  are  very  sensitive  to  the  high  frequencies  of  the  gravity  field  and  thus  should  be 
computed  with  the  densest  grid  of  data  available.  Presently  dense  grids  of  gravity  data  required 
for  computing  the  correction  terms  accurately  are  available  only  in  die  United  States  and  in  some 
countries  in  Westom  Europe.  Therefore  consideration  of  conection  terms  beyond  the  gi  term  in 
the  computation  of  the  height  anomaly  all  over  the  globe  seems  impossible  as  of  today  (1994). 

Non-consideration  of  the  g2  term  in  the  hei^t  anomaly  may  lead  to  an  error  of  about  2 
cm  in  a  moderately  undulating  terrain  and  more  in  rough  mountainous  areas.  Sideris  and 
Schwarz  (1986)  have  pointed  out  diat  the  use  of  free-air  gravity  anomalies  corrected  for  terrain 
(also  called  as  Faye  anomalies)  in  the  computation  of  the  height  anomaly  makes  the  gn  terms 
smaller,  smoother  and  easier  to  interpolate.  Hence  the  use  of  Faye  anomalies  in  the  computation 
can  minimize  the  error  due  to  neglecting  the  g2  term. 

Now  considering  only  up  to  Molodensky's  cmrecdon  term  gi,  equation  (3.25)  simply 
becomes: 


;,=C.+^JKAgS+*.)«'t')<K»  (3-30) 

If  terrain  corrected  free-air  anomalies  are  available,  then  from  equation  (24)  of  Wang  (1993),  we 
can  rewrite  equation  (3.30)  as: 

Cp=Co+^JJ(AgS+C)S(V)(to-^H;-:^H?-:^8h>  (3.31) 

where  ( Ag^ + C)  is  the  gravity  anonudy  cmiected  for  terrain  and  other  systematic  effects  listed 
in  Section  3.1.1, 


Hp  is  the  normal  height  of  the  point  P  at  which  the  height  anomaly  is  computed. 


5h^  =  0.453  •  0.018  sin<t>  +  0.087  cos<)>  cosX  -i-  0.204  cos<^  sinX,  (km^)  where  0,  X.  are  the 
geographic  latitude  and  longitude  of  the  point  P,  respectively,  and 

Aga  =(^gp  -27iGpHp)  is  the  Bouguer  anomaly,  Agp  is  the  surface  free-air  anomaly,  G  is 
Newtonian  gravitational  constant  and  p  is  the  density  of  topographic  mass. 

Equation  (3.31)  is  based  on  the  assumption  that  the  gravity  anomaly  Ag^  is  linearly  correlated 
with  elevation,  and  it  provides  an  accurate  solution  for  computing  the  height  anomaly  when 
terrain-corrected  free-air  anomalies  are  available. 

Combining  equations  (3.23)  and  (3.31),  we  can  write: 

Y  Y  Y 


iq'CQio//S(\|/)daj 


(3.32) 


where 


In  equation  (3.32)  n  denotes  the  number  of  regional  vertical  datums  considered  in  the  Global 
Vertical  Datum  definition. 

From  Figure  3.2,  we  note  that  approximately. 


(3.33) 


where  Hp  is  the  normal  height  of  the  station  referred  to  the  global  vertical  datum. 


H*  is  the  normal  height  of  the  station  referred  to  the  regional  vertical  datum,  and 

^=w„-wi 


also 

hp=H*p+Cp  (3.34) 

where  hp  is  the  ellipsoidal  height  of  the  point  P. 

Substituting  =  — p  —  ^  °  in  equation  (3.32)  where  5(GM)  =  GM  -  (GM)o  is  the 

KY  Y 

difference  between  the  true  geocentric  gravitational  constant  and  that  of  adopted  reference 
ellipsoid,  we  get: 
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(3.35) 


C,  =  4I7 1 JJK +C)'s(v)do,  tq'C,,.J/S(VXi«», 

_.^h;  -^2£hp  -55£5h* 

y  y  y 

Since  the  GM  value  is  known  to  an  accuracy  of  mproximately  0.003  ppm,  GM  shall  be  regarded 
as  known  in  the  discussion  and  one  can  assume  ^GM)  =  0  (^pp  and  Balasubramania,  lw2)  in 
equation  (3.35). 


From  equation  (3.34)  we  can  write: 

^p  =  hp-Hp 

Comparing  equations  (3.35)  and  (3.36): 
.  „._-AW„  .  R  “ 

lip  “  rtp  “ 


H;=-^+^l/J(Ag^+cyS(v)da,-^IqiCQ^op(¥)da3 


-.^h;  — Z^Hp  -i59£5h* 
Y  Y  Y 


From  equation  (3.33)  H*  =  H*  + 


CqjO 


Substituting  (3.33)  in  (3.37)  we  get: 


Hr=^+i^+5^t{f(AgS+C)'S(v)dOi-;iiq‘C,„|fS(v)d. 


_^h;  -^^Hp  -i^6h2 
Y  Y  Y 


Substituting  (3.24')  in  (3.38)  and  rewriting  the  known  and  unknown  terms,  we  get: 

h,  -h;  |/j(Ag£  +c)'  s(v)dOi +^h; +i^Sh» 

ijJS(y)dOi 


AW„  .  Cq,.  ,  2-.^ 


(3.36) 


(3.37) 


(3.38) 


(3.39) 


In  equation  (3.39)  the  left  hand  side  can  be  computed  using  the  observaticms.  Let  it  be  (tenoted 
byY.  Then, 


Y  =  hp-H*- 


R 

4jcy 


i//(Ag£  +cys(y)dOj+^H;+:^H?+5^8h’  (3.40) 

'  '  T  T  T 


where  hp  is  the  ellipsoidal  height  of  the  point  P,  and  H*  is  the  normal  height  of  the  station 
referred  to  the  i^  vertical  datum  in  which  it  lies. 
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Since  in  equation  (3.40)  the  terms 


Aga  nGp 

Y  *  y  ’ 


etc.  are  small,  and  the  difference  between 


H*  and  Hp  is  also  very  small,  we  can  replace  Hp.  Hp  in  the  right  hand  side  of  equation  (3.40) 
by  H*  and  H^,  respectively. 


Denoting  the  gravimetric  height  anmnaly  (without  zero  height  anmnaly  term)  by  Cl.  we  can  write 
equation  (3.39)  as: 


(3.41) 


(3.42) 


Equation  (3.41)  is  the  required  observation  equation  when  free-air  anomalies  in  Molodensky's 
sense  and  normal  heights  of  the  station  above  regional  vertical  datum  are  given. 

3.2.2  Procedure  When  Free-air  Anomalies  in  aassical  Sense  and  Orthometric  Heights  of  the 
Stations  are  Known 

From  equation  (2.S)  we  know  that  the  orthometric  height  of  a  surface  point  P  with  respect 
to  an  equipotential  surface  W  «  Wq  (in  our  case  the  global  vertical  datum)  is  given  by: 

g 

where  Cp  is  the  geopotential  number  at  point  P  (Cp  =  Wq  -  Wp)  and  g  the  nsan  gravity  along 
the  curved  vertical  between  the  surface  point  and  its  corresponding  point  Pq  on  the  global 
vertical  datum  (Refer  to  Figure  3.3).  In  the  same  way  the  ortiiometric  height  of  the  point  P  with 
respect  to  the  regional  vertical  datum  'i'  (W  =  WO  can  be  written  as: 

g 


where  Cp,  =  W{,  -  Wp 


(3.43) 


Geoptilmmgh  point  P  Ws^ 


W  =  Wo 
W  =  Wo 


U  =  Uo 


Figure  3.3  Various  Equipotential  Surfaces  and  Orthometric  Height 

In  Figure  3.3,  Qi  can  be  considered  as  the  fundamental  benchmark  of  the  i-th  regional 
vertical  datum.  The  potential  of  the  reference  surface  for  the  global  vertical  datum  is  Wq.  We 
define  as  before  Cq.o  =  Wo  -W^  as  the  potential  difference  between  the  equipotential  surfaces  of 
the  regional  datum  i  and  the  global  vertical  datum.  Introducing  an  arbitrary  ellipsoid  where  Up  is 
the  potential  on  the  surface,  we  define  AWq  =  Wq  -  Uq.  Also  in  the  above  figure  the  distance 
between  the  global  vertical  datum  and  the  arbitrary  ellipsoid  is  given.  This  depends  on  the 
disturbing  potential  at  the  point  Po.Tp^ ,  and  AWq  =  Wq  >  Uq. 

From  equations  (3.1)  and  (3.4)  we  know  the  free-air  anomaly  in  a  classical  approach  including 
terms  only  up  to  the  second  order  is  written  as: 


Agc  = 


-Ta 


Qo 


(3.44) 


where  H  is  the  orthometric  height  and  is  the  normal  gravity  on  the  ellipsoid.  Since  the 
oithometric  heights  of  the  stations  are  knowm  with  res^t  to  the  regional  vertical  datum  i,  H>,  the 
gravity  anomalies  computed  using  equation  (3.44)  will  refer  to  the  regional  vertical  datum  only. 
Denoting  the  free-air  anomalies  referring  to  the  regional  vertical  datum  i  as  Agj^,  we  can  derive 
the  relationship  between  Ag^  and  Ag^  that  refer  to  the  global  vertical  datum  as  follows: 

For  an  arbitrary  reference  ellipsoid,  following  a  simple  generalization  of  Brun's  formula 
we  can  write  (refer  to  Figure  3.3) 


N 


(3.45) 


where  Tp^  refers  to  the  disturbing  potential  at  the  point  Pq  on  the  global  vertical  datum. 
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Similarly  for  the  regional  vertical  datum'!'  we  can  write 


(3.46) 


where  Tn  refers  to  the  disturbing  potential  at  the  point  on  the  regional  vertical  datum,  and 
AWj,  s  -  Uo  is  die  potential  difference. 

Frran  (3.45)  and  (3.46)  we  get 


N -N‘ -  7[K -T,,) -(AW,  -  AWi)] 


since  »  T^,  we  can  write 

y 

with  eqn.  (3.47),  eqn.  (3.45)  can  be  modified  as 

y  y 

From  equation  (2-148')  in  Heiskanen  and  Moritz  (1967)  we  have  the  relation: 

ar. 


Ag.  =  — gi+^N 
cm  cm 


(3.47) 


(3.48) 


(3.49) 


where  Age  is  the  gravity  anomalies  defined  in  a  classical  sense  referring  to  the  global  vertical 
datum, 

N  is  the  geop-spherop  separatitm  between  the  GVD  and  the  arbitrary  reference  ellipsoid, 
y  is  the  normal  gravity,  and 

~  is  the  derivative  along  the  ellipsoidal  normal, 
on 

Similarly  for  the  case  of  the  regional  vertical  datum '!',  we  can  write 


N‘ 


Now  substituting  for  N‘  from  equation  (3.48),  we  get 


(3.50) 


(3.51) 


In  spherical  approximation,  setting  r  =  R, 
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(3.52) 


Substituting  (3.52)  in  (3.51) 

From  eqn.  (2-179')  in  Heiskanen  and  Moritz 

From  equations  (3.53)  and  (3.54)  we  get 


Age  —  Agi  +"^0^0 


(3.53) 


(3.54) 


(3.55) 


From  equation  (2-181)  in  Heiskanen  and  Moritz  (1967),  we  can  write: 

N(P)=5^-^+J^|AfcS(¥)d<j  (3J6) 

Substituting  (3.47)  and  (3.55)  in  (3.56)  and  assuming  that  5(GM)  is  equal  to  zero  because  of  the 
high  accuracy  in  the  determination  of  GM  today  (1994)  as  noted  earlier,  we  get: 

N‘<P) = + |c,,,]dOj  (3.57) 

where  'i'  refers  to  the  i<^  regional  vertical  datum  and 

'j'  =  1, 2...n,  n  being  the  number  of  regional  vertical  datums  ctmsidered. 

The  free-air  anomaly  used  in  equation  (3.57)  should  be  conected  for  the  systematic  errors 
listed  in  Section  3.1.1  and  also  for  the  terrain  effect  If  Faye  anomalies  are  used  in  the  Stems' 
formula  in  equation  (3.57)  then  by  Helmen's  second  method  of  condensation,  an  indirect  effect 
correction  is  needed  to  be  applied  to  the  computed  undulation.  The  required  equation  for 
computing  the  indirect  effect  as  given  by  Grushinsky's  formula  (Wichiencharoen,  1982)  is: 

5N'  =  (independent  of  j  =  l,...n)  (3.58) 

where  the  quantities  G,  p,  y  ate  defined  eariier  and  H>  is  the  height  of  the  station  in  the  regional 
vertical  datum. 

Assuming  the  fiee-air  anomalies  Age  are  corrected  for  systematic  corrections  and  also  for  terrain, 
equation  (3.57)  gets  nxxiified  as: 
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N'(P>=^+^+^tj|S(v)[(Agi+cy+|c,,,]<lo,-i^H‘‘  (3.59) 

where  (Ag^  +  C)^  refers  to  the  fiee-air  anomalies  defined  in  a  classical  sense,  conected  for  terrain 
and  other  systematic  corrections  lefeiring  to  regional  vertical  datum  'j',  for  j  =  1, 2,  ...n  wtere  'n' 
is  the  number  of  regional  vertical  datum  considoed  in  the  global  vertical  datum  definition. 

Assuming  that  we  have  a  set  of  space  geodetic  stations  on  regional  vertical  datum  'i',  for 
which  we  know  the  ellipsoidal  heights  li'  above  an  arbitrary  geocentric  reference  ellipsoid  and 
the  orthometric  heights  H‘  above  the  local  vertical  datum,  then  we  can  write  from  equation  (2.1) 
and  Hg.  3.3: 

N*  =  h  -  H»  (3.60) 


comparing  equations  (3.59)  and  (3.60) 
h-H‘  = 


-AW„  C 


o+i:Qii+ 


3|^ijJS(v((Agj4C)Vlc„.]*,,-5feH‘’ 


t|S(V)(Ag'  +C)‘dai  +2  ic 


TCGp 


Hf 


Re-arranging  the  known/observable  and  unknown  terms,  we  get: 


>■  -  H'  iJJ  (Ag:,  ■<-c)'s«i(Xto, 


Denoting  the  left  hand  side  of  equation  (3.61)  as  Y,  we  write: 


(3.61) 


Y  =  h-H‘ 


(3.62) 


The  third  term  on  the  right  hand  side  of  the  equatimi  refers  to  the  geoid  undulation  computed 
using  terrain  corrected  gravity  anomalies  without  the  zero  undulation  (Nq)  term.  Denoting  the 
third  term  by  Nj  and  from  equation  (3.61)  and  (3.62)  we  can  write: 


(3.63) 


where  Ni  =  ;^||(Ag-+C)'s(v)<tej-^Hf 


(3.64) 
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Equation  (3.63)  is  the  required  observation  equttion  when  free-air  anomalies  in  the  classical 
sense  and  orthometric  heists  of  the  station  aboW  the  regional  vertical  datums  are  given. 

3.2.3  Procedure  When  Free-air  Anomalies  in  aassical  Sense  and  Normal  Orthonaetric  Heights 

of  Stations  arc  Known 

As  explained  in  section  2. 1.2.2,  nmmal  mhometric  heights  are  used  in  countries  which 
lack  dense  gravity  data  coverage.  However,  they  are  consider^  for  all  practical  purposes  as  a 
coarser  approximation  to  true  o^ometric  heights.  Therefore,  the  procedure  described  in  Section 
3.2.2  can  be  used  for  setting  up  the  requir^  observation  equations  in  diis  case  also.  Since 
normal  onhometric  heights  are  not  true  orthometric  heights,  such  heights  can  produce  a 
systematic  error  in  the  definition  of  global  vertical  datum. 

3.3  Computation  of  Gravimetric  Height  Anomalv/Unduladon 

To  achieve  the  highest  accuracy  in  the  coinputation  of  gravimetric  height  anomaly,  Ci  •  in 
equation  (3.42)  or  the  gravimetric  undulation,  Nj,  in  equation  (3.64),  the  pr«%dure  that  can  be 
followed  is  to  combine  the  surface  gravity  antnnaly  data  defined  in  a  regional  vertical  datum 
limited  to  a  small  cap  (100  to  200  km)  about  the  computation  point  and  coefficients  from  an 
accurate  geopotentiai  iiKxlel.  Some  of  the  currently  used  techniques  for  computing  gravimetric 
height  anomaly/geoid  undulation  are  discussed  in  tlw  following  p^graphs. 

3.3. 1  Modified  Stokes'  Technique 

This  technique  makes  use  of  gravity  anomaly  data  in  a  cap  around  the  point  of  interest 
and  a  set  of  potential  coefficients  from  an  accurate  global  geopotentiai  model  to  compute  the 
height  anomaly/geoid  undulation.  This  technique  has  been  discussed  by  several  authors;  to 
mention  a  few,  Despotakis  (1987),  Pavlis  (1991),  and  Sjdberg  (1991).  From  equation  (2.40)  in 
Despotakis  (1987),  we  can  write  die  modifled  Shoes’  equaticni  (here  incorporating  the  indict 
effect  term): 


N=3|^JjAgS(yXlo+^_|Q.(y)ag,(e.>.)-5^H>  (3.65) 

where  Ag  are  the  free-air  gravity  anomalies  corrected  for  systematic  errors  listed  in  Section  3.1.1 
and  considered  in  a  cap  radius  yc> 

Qn  are  the  Molodensky  truncation  coefficients  defined  by  equation  (7-34)  in  Heiskanen 
and  Moritz  (1967),  and 

Agn  (6,  X)  the  n<b  degree  harmonic  representation  of  Ag  (6,  X). 

The  frrst  term  on  the  tight  hand  side  of  equation  (3.65)  corresponds  to  the  contribution  of  gravity 
anomalies,  in  a  cap  (y  s  y^),  to  the  undulaticm  computation  and  the  second  term  to  the  remote 
zone  contribution  to  the  un^ation  computation.  When  terrain  corrected  gravity  anomalies  are 
used,  equation  (3.65)  gets  nxxiified  as: 

R  R^  -  Irf-n 

oc  zya-o  y 
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where  AgJ (6.X)  is  the  n*!*  degree  harmonic  term  of  (Ag  C). 

From  discussions  in  Raj^  (1993a),  Ag][(6,X)  in  the  above  equation  is  related  to  Agn  (0,  X)  in 
equation  (3.6S)  as  follows: 

AgI(e.X)  =  Ag,(0.X)+C.(0.X)  (3.67) 

where  C  is  the  terrain  correction,  the  coefficients  of  C  can  be  found  by  replacing  Ag  (0,  X)  in 
equation  (4)  (Rapp,  1993a)  with  C,  the  terrain  correction. 

Since  potential  coefficient  models  available  for  use  have  yielded  estimates  of  the 
potential,  we  can  calculate  AgI(0,X)  from  potential  coefficients  as: 

AgI(0,X)  =  ^!^(n-l)^^j  i^(C„cosmX+S„sinmX)P„(cos0)  (3.68) 

Using  (3.66)  and  assuming  gravity  anomalies  used  in  the  ctq)  refer  to  only  one  datum  'i'  in  which 
the  station  lies,  from  equation  (3.42)  we  get: 


Ci  =  ^//(AgS  +Cj  S(v)do.  +  A|Q.(y)AgJ(e,X)- 
i  Y  ‘  Y  ‘  Y  J 


(3.69) 


Similarly,  using  (3.66)  and  assuming  gravity  anomalies  used  in  the  cap  refer  only  to  the  datum  'i' 
in  which  the  station  lies,  ffimi  equation  (3.64)  we  get: 


Ni  =  ^Jj(AgJ  +cJs(y)do,+AiQ.(v)Ag:(eA)-J^H?  (3.70) 

The  gravimetric  height  anomaly  and  the  gravimetric  undulation  Nj  computed  by 
equations  (3.69)  and  (3.70)  with  the  Modified  Stokes'  technique  are  based  on  sphericsd 
approximation.  But  to  get  the  solution  with  respect  to  an  ellipsoid  as  reference  surface  we  need 
to  apply  ellipsoidal  corrections  to  the  height  anmnaly/undulation  computed.  Rapp  (1981)  has 
develof^  procedures  to  compute  the  correction  need^  to  fully  refer  the  solution  to  an  ellipsoid 
surface  (equation  (31),  ibid.).  Despotakis  (1987,  C2iapter  VI)  has  also  derived  the  ellipsoid 
correction  term  to  be  applied  to  the  height  anomaly/undulation  computed  using  gravity 
information  in  a  cap  surrounding  the  computatitHi  point  combined  with  potential  coefficient 
information,  following  the  similar  lines  as  in  Rapp  (1981).  Denoting  the  ellipsoidal  correction 
term  as  AN,  equations  (3.69)  and  (3.70)  get  modif^  as: 

Q  =^//(Ag£  +cjs(v)<to. iQ.(y)AgI(e,)i)- 

.  (3.71) 

+^5h*  1+AN 
I  Y  Y  Y  J 

Ni  Jj(Ag?  +cJS(v)da.  +  A  £Q,(y)AgI(e,X)-5^Hf  +  AN  (3.72) 
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(3.73) 


AN=:|i  i(Q„ -X.)  t  (G„cosmX+H«sinmX)P«(sin^)+c47-|sin*4>]No 

■*T  B»0  m^O  L4  4  J 

where  are  Molodensky's  truncation  coefficients  (Heiskairen  and  Moritz  (1967,  Section  7.4)), 

Xn  are  the  Fourier  coefficients  of  Stokes'  function  beyond  the  cap  radius  Vc>  and 

Gnin>  Hfim  are  the  coefficients  defined  in  equation  (39-81)  in  Moritz  (1980). 

The  first  term  in  equation  (3.73)  arises  from  the  fact  that  the  ellipsoidal  effect  e^Agl* 
where  Agl  is  defined  by  equation  (39-80)  in  Moritz  (1980),  has  to  be  remov^  from  the  terrestrial 
gravity  anomalies  Age  us^  in  the  cap.  As  shown  in  Despotakis  (1987)  this  correction  term  is 
dependent  on  the  ctq)  radius  yc  in  which  the  gravity  anoiiudy  data  are  used  and  its  magnitude 
amounts  to  about  1-2  cm  when  yc  =  2**. 

The  second  term  in  equation  (3.73)  corresponds  to  the  ellipsoidal  correction  to  the 
spherical  undulation  itself  and  is  given  as  equadon  (39-21)  in  Moritz  (1980)  also,  litis  term  has 
b^n  derived  assuming  average  normal  gravity  being  used  at  the  station.  If  actual  normal  gravity 
computed  at  the  stations  are  used  in  computing  the  undulation,  this  term  becomes  zero. 

3.3.2  Least  Squares  Collocation  Technique 

As  discussed  in  section  16  of  Moritz  (1980),  Least  Squares  Collocation  (LSC)  technique 
provides  the  most  accurate  result  obtainable  on  the  basis  of  ^ven  data  when  pritMr  information  is 
included.  Prediction  of  hei^t  anomaly  Cl  or  gravimetric  undulation  N| ,  based  on  the  corrected 
free-air  anomalies  given  in  a  regional  datum  can  be  done  accurately  using  least  squares 
collocation  technique  provided  the  requited  covariances  are  known. 

The  statistical  approach  to  collocation  is  based  on  the  theory  of  stochastic  processes  in 
which  the  gravity  anomaly  field  is  treated  as  a  stadonary  stochasdc  process  on  a  sphere  with  the 
assumpdon  that  both  the  observations  and  signals  represent  ergodic  processes  with  zero 
expectation  everywhere  (Moritz,  1980,  p.  76)  and  the  average  product  of  two  signal  quantities  is 
the  covariance  between  them.  Therefore,  in  the  case  of  gravity  anonudies  it  is  assumed  that, 

M{Ag)  =  0  =  E{Ag)  and  M{Ag(P)  •  Ag(Q))  »  C  (P.  Q)  =  E{ Ag(P)  •  Ag(Q)}  (3.74) 

The  definition  of  the  mean  operator  M('}  as  a  honoogeneous  and  isotropic  average  over  the 
whole  sphere  is  an  important  assumption  in  this  technique. 

On  a  global  scale,  the  anomalous  gravity  field  exhibits  a  seemingly  random  behavior  necessary 
for  a  statistical  interpretation  of  the  freld.  But  wten  the  area  of  interest  is  of  limited  extent,  as  in 
our  case  with  limited  surface  gravity  data,  the  global  gravity  field  behavior  represents  a 
systematic  component  and  short  wavelength  features  become  prominent  from  the  statistical  point 
of  view  (refer  to  Figure  3.4  given  as  Figures  1(a)  and  1(b)  in  Collier  and  Leahy,  1992). 
Therefore,  the  conditions  for  equation  (3.74)  need  to  be  modified  with  the  requirement  that  the 
local  average  behavior  must  be  zero  and  a  local  covariance  function  must  be  us^. 


Hgure  3.4  Gravity  Field  <m  a  Global  and  Local  Scale 
(from  Q)llier  and  Leahy,  1992) 

To  satisfy  the  requirement  of  zero  expectation  as  well  as  that  of  isotropy  and 
homogeneity  for  the  covariaiice  function,  it  is  essential  that  the  presence  of  all  systematic  trends 
in  the  local  behavior  of  the  gravity  field  be  eliminated  from  the  gravity  data  prior  to  generating 
an  empirical  covariance  function.  After  the  removal  of  the  trend  component,  the  remainder  of 
the  gravity  field  shall  here  be  considered  to  behave  in  a  stochastic  fashion.  Global  geopotential 
n^els  up  to  degree  and  order  360  can  be  used  for  removing  the  trend  in  the  gravity  field. 
Practical  computation  of  reference  gravity  ammialies  from  a  set  of  potential  coefficients  (up  to 
degree  360),  computation  of  a  globid  covariance  function  and  then  scaling  it  to  fit  local  data  are 
discussed  in  detail  in  BaSi6  and  Rapp  (1992). 

From  equation  (18-21)  in  Moritz  (1980),  the  required  expression  for  predicting  the 
residual  geoidal  undulation/height  anomaly  at  the  space  geodetic  stations,  using  the  LSC 
technique  is  given  by 

s  =  C„C-‘Y  (3.75) 

where  s  is  the  predicted  residual  geoid  undulation/height  anomaly, 

Y  are  the  residual  gravity  anomalies  obtained  sbher  removing  the  reference 
anomalies  computed  using  global  geopotential  model  from  the  observed  free-air 
anomalies  conectod  for  terrain  and  ot^  systematic  effects, 

Cst  is  the  cross  covariance  matrix  between  unobserved  and  obsarved  quantities,  and 
"C  is  sum  of  the  covariance  matrices  of  the  signal  t  and  noise  n  (Cu+Ciiii>=D{  Y). 

In  addition  to  providing  an  accurate  prediction  of  height  anomaly/geoid  undulation,  least 
squares  collocation  technique  can  also  provide  the  accuracy  of  the  predict^  quantities  which  can 
be  computed  using  the  error  covariance  matrices  given  as  equation  (17-44),  (17-47)  and  (17-48) 
in  Mmitz  (1980). 
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E-*C,-C.C-‘C.  (3.76) 

where  Eat  is  the  enor  covariance  matrix  of  {mdicted  geoid  undulatimi, 

Cgs  is  the  auto  covariance  matrix  of  unobserved  gec^  undulation. 

Cat.  Cta  are  the  cross  covariance  matrices  as  before  Cat=C(s,Y}.  and 

C  is  the  total  covariance  matrix  of  observations  (Cu  +  C^sD{  Y). 

The  height  anomaly/geoid  undulaticms  predicted  using  LSC  technique  are  to  be  cmrected 
for  ellipsoidal  corrections  following  the  procedure  described  in  Section  39  of  Moritz  (1980).  The 
discussions  given  in  there  indicate  that  tte  ellipsoidal  conection  conq}uied  and  applied  in  Section 
3.3.1  for  the  Modified  Stokes'  Technique  can  be  applied  in  this  case  also.  As  pointed  out  earlier, 
the  magnitude  of  ellipsoidal  correction  is  in  the  ot^  of  1*2  cm.  The  geoid  u^ulation  predicted 
using  LSC  technique  should  also  be  corrected  fcH*  indirect  effect  discussed  in  Section  3.2.2  (eq. 
(3.58)). 


The  practical  problem  in  using  this  procedure  (mainly  with  available  software)  may  be 
the  inversion  of  the  matrix  C  in  equation  (3.75)  as  well  as  the  equadtm  (3.76).  The  size  of  the 
matrix  to  be  invened  will  be  (nxn)  where  n  is  the  number  of  observations.  Procedures  to 
circumvent  the  problem  of  inverting  the  matrix  have  also  been  developed  for  predicting  using  the 
LSC  technique.  However,  depending  on  the  limitations  of  the  computers  us^  in  terms  of  CPU 
time  and  memory,  the  number  of  observations  used  for  predicting  will  be  curtailed. 

3.3.3  Fast  Fourier  Transform  (FFDTechniQue 

The  conversion  of  gravity  anomalies,  defined  in  a  regional  vertical  datum  to  geoid 
heights  can  be  done  using  FFT  method  also,  after  removing  the  long  wavelength  part  of  the 
gravity  field  using  an  accurate  geojwtential  model,  complete  to  degree  and  order  360.  This 
method  takes  advantage  of  the  possibility  of  expressing  the  planar  approximation  of  Stokes' 
integral  in  a  two  dimensional  convolution  form.  According  to  the  convolution  theorem 
(Bracewell,  1978),  the  Fourier  transform  of  a  convolution  is  the  poduct  of  the  spectra  of  the  two 
functions  being  convolved.  Detailed  discussions  on  expressing  the  Stokes'  integral  in  a  planar 
form  and  its  transfmmation  to  sptxtral  form  are  given  in  Sections  2.3  and  2.5  of  Zhao  (1989). 
The  planar  s^proximation  underlying  this  approach  is  peimissible  due  to  the  use  of  a  high-order 
spherical  harmonic  reference  field,  yielding  near-negligible  indirect  effects  (ForsbCTg,  1990). 
However,  it  is  also  possible  to  apply  HFT  directly  to  the  spherical  Stokes'  formula.  This  gives 
more  accurate  results  and  therefore  the  extensiwi  of  the  area  of  integration  can  be  enlarged 
(Strang  van  Hees,  1990).  A  review  of  FFT  methodology  may  be  found  in  Schwarz  et  al.  (1990) 
and  the  different  sources  of  error  in  die  use  of  FFT  mediods  for  conversion  of  gravity  anomalies 
to  geoid  heights  in  Farelly  (1991). 

Though  the  high  computation  speed  of  FFT  makes  this  technique  a  very  efficient  (me  in 
various  scientiric  field^,  computation  of  height  anomalyAmdulation  at  one  station  (not  falling  at 
any  grid  coma-  of  given  gravity  anomalies)  does  not  give  any  speciHc  advantage  for  this  method 
Accuracy  of  prediction  is  almost  the  same  as  modified  Stcdtes'  technique  or  even  less  due  to 
interpolation  of  the  final  geoid  undulation  to  the  station  from  grid  comers,  where  N  are  predicted 

Most  of  the  regional  geoid  height  models  such  as  GEOID93  (Milbert,  1993),  GSD91 
(Mainville  and  V6ronneau,  1990),  and  I*^G-89  (Fcnsberg,  1990)  makes  use  of  FIT  technique  for 
their  development.  The  required  values  for  Cl  and  Nj  for  ^nations  (3.42)  and  (3.64)  can  also  be 
interpolated  firom  the  regional  geoid  height  tnodels  where  it  is  available. 
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The  discussions  in  this  section  follow  the  procedure  suggested  by  Rapp  and 
Balasubramania  (1992). 

From  equations  (3.41)  and  (3.63),  we  get  the  observation  equations  at  each  of  the  space 
geodetic  stations  'k'  where  k  =  1,  2...z  and  ‘z*  is  the  total  number  of  stations  used  in  the  global 
vertical  datum  deHnition: 


-AWo  .  C 


(3.77) 


we  now  define  the  integral: 


W  =  ijJS(V)d<Tj 


(3.78) 


Considering  terrestrial  gravity  anomalies  only  around  a  small  cap  (\|f0  =  2**)  we  write: 
J(Vc)  =  4Ts(V)sinvdv 


(3.79) 


Equation  (3.79)  is  also  equation  (2-225)  in  Heiskanen  and  Moritz  (1967).  The  integration  is 
given  by  Lambm  and  Darling  (1936,  p.  103).  We  have: 

J(Vc)  =  jl +4sin'^  -cos\|fj  -6sin^  — ^sin*  (3.80) 

Now  assuming  that  the  anomalies  in  the  cap  refer  to  only  one  vertical  datum  as  discussed 
in  Section  3.3.1,  in  which  the  station  lies,  then  for  j  =  i  equation  (3.77)  reduces  to: 


Yk=-^+^(l+2J(vJ) 


(3.81) 


If  gravity  anomalies  in  other  vertical  datums  are  involved,  additional  terms  are  needed  in 
equation  (3.81). 


The  general  form  of  equation  (3.81)  can  be  written  as: 
Y=AX  +  e 


(3.82) 


Where  Y  is  the  misclosure  vector  from  equations  (3.41)  or  (3.63),  A  is  the  coefficient  matrix,  X 
are  the  parameters  of  the  model  (AWo  and  CQjo),  and  e  the  observation  error  vector,  l^e 
coefficient  of  AWq  will  be  -1/y  and  the  coefficient  of  Cqjo  be  (l+2J(v))/y,  where  we  assume 
the  same  cap  size  for  all  stations. 

In  our  global  vertical  datum  problem,  there  will  be  z  observation  equations  and  I+l 
parameters  where  I  is  the  number  of  vertical  datums  considered.  In  this  system  there  will  be  a 


rank  (tefect  of  1  because  nothing  has  been  set  to  define  a  fixed  datum  or  surface.  This  defect  can 
be  removed  by  introducing  a  constraint.  One  such  procedure  as  discussed  by  Rapp  and 
Balasubramania  (1992)  is  to  carry  out  a  free  network  adjustment  as  is  done  for  some  ty^s  of 
horizontal  netwoi^s  or  to  introduce  a  physically  meaningful  constraint  on  the  system.  We  define 
the  constraint  equation  such  that: 


X(Wo-Wi,)  =  I(Cq,o)  =aminimum 


iBl  ial 

To  implement  equation  (3.83)  we  introduce  a  datum  constraint  such  that 
KX  =  0  with  K:=[0. 1, 1 1] 

where  K  is  the  Ixm  coefficient  matrix  of  the  constraints  with  rank  K  =  1,  and 

*A1 


(3.83) 


(3.84) 


rank 


=  m  =  I  +  1. 
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Defining  the  Lagrange  target  function  as: 

4>=(Y-AX)TI;'(Y-AX)+2X(KX)  (3.85) 

a  least  squares  solution  based  C'  equation  (3.85)  leads  to  the  following  system  of  equations: 

(3.86) 

The  solution  of  the  above  equation  yields  the  following  estimated  parameter  vector 

x=(a^i;'a+ktk)'‘at5:;‘y  (3.8?) 

The  error  covariance  matrix  of  the  estimated  parameters  would  be: 

li  =  (A^  I;‘  A + K^K)’'  A^  I;‘  A(  A’'  I;‘  A + (3.88) 


The  elements  of  K  are  in  principle  arbitr^.  The  simplest  selection  of  the  elements 
would  be  to  set  all  but  one  element  to  1  as  shown  in  equation  (3.84).  This  would  imply  that  all 
vertical  datums  would  contribute  equally  to  the  determination  of  the  reference  surface  of  the 
global  venical  datum.  An  alternative  solution  would  be  one  where  some  vertical  datums  would 
be  given  greater  weight  in  defining  the  global  vertical  than  others.  Proper  solution  of  the  actual 
numerical  values  fOT  the  K  elements  depends  on  various  factors  such  as  the  number  and  accuracy 
of  space  geodetic  stations  considered  in  each  datum,  number,  type  and  accuracy  of  tide  gauges 
consider^  in  (fefining  the  regional  vertical  datum  etc.  Also  we  note  that  Wq  does  not  (hrectly 
enter  the  constraint  defined  by  equation  (3.84).  Therefore,  the  K  value  corresponding  to  AWq 
parameter  must  be  chosen  as  zero. 

Equations  (3.87)  and  (3.88)  are  the  key  equations  in  estimating  the  adjusted  parameters 
defining  the  global  vmical  datum  and  their  accuracy  estimates.  Actual  estimation  of  parameters 
and  their  accuracy  are  carried  out  in  the  next  chapter  based  on  available  data  in  various  regional 
vertical  datums. 


CHAPTER4 

REALIZATION  OF  THE  FIRST  ITERATION  GLOBAL  VERTICAL  DATUM 


Based  on  the  modeling  procedures  discussed  in  Chapter  3,  a  Hrst  iteration  test 
computation  to  estimate  the  parameters  defining  the  Global  Vertical  Datum  was  carried  out  with 
available  data  from  six  regional  venical  datums.  The  different  regional  vertical  datums 
considered  for  realizing  the  Global  Venical  Datum  (GVD)  are:  North  American  Vertical  Datum 
88  (NAVD88),  Ordnance  Datum  Newlyn  (England),  IGN  69  (France),  NN  (Germany),  AHD71 
(Australia),  and  Scandinavian  Datum.  Data  used  in  the  numerical  investigations,  corrections 
applittl  to  them,  different  techniques  used  in  evaluating  certain  gravimetric  quantities  and  the 
results  obtained  are  discussed  in  the  following  sections. 

4J  Data.Usfid  in  Diffgrept  Regional  Vertical  Datums 

The  ITRF91  coordinate  system  was  selected  as  the  consistent  coordinate  system  for 
defining  the  station  cotndinates  of  the  space  geodetic  stations  such  as  VLBI,  SLR,  GPS,  etc., 
used  in  the  realization  of  GVD.  Station  cooi^nates  in  the  ITRF91  (for  epoch  1988.0)  for  all 
space  stations,  except  for  one,  used  ir  this  study  were  obtained  firom  Table  T6  in  lERS  Technical 
Note  12  (Boucher  et  al.,  1992(a),  pp.  81 -86).  An  ideal  ellipsoid  with  the  following  parameters 
was  adopted  as  the  reference  ellipsoid  for  this  study: 

a  =  6378136,3  m 
f  =  1/298.257222101 

GM  =  3986005  x  108m3s-2  (4- 1 ) 

m  =  72921 15  X  10-11  rad  s-1 

The  last  three  parameters  are  the  same  as  the  parameters  defining  the  GRS80  ellipsoid.  The  list 
of  stations  us^,  their  ellipsoidal  coordinates,  ellipsoidal  heights  and  heights  in  local  vertical 
datum  are  given  in  Table  4.1. 

4.1.1  North  American  Vertical  Datum  88  (NAVD88) 

This  vertical  reference  datum  has  been  recently  established  (Zilkoski  et  al.,  1992)  based 
on  a  irnnimum-constraint  adjustment  of  Canadian-Mexican-US  leveling  observations,  holding 
fixed  the  height  of  the  primary  tidal  bench  mark,  referenced  to  the  locid  mean  sea  level  height 
value  at  Father  Point/Rimouski,  Quebec,  Canada.  NAVD88  values  are  given  in  Helmert 
orthometric  height  units.  This  height  system  will  eventually  replace  the  National  Geodetic 
Vertical  Datum  29  (NGVD29)  height  system  used  in  the  United  States  during  the  next  five  to 
seven  years.  The  NGVD29  height  system  was  established  by  constraining  the  local  mean  sea 
level  at  26  tide  gauge  stations  (21  in  the  US  and  5  in  d^ada)  held  fixed  at  0.0  elevation.  Die  list 
of  distortions  present  in  the  NGVD29  height  system  are  given  in  Table  1  of  Zilkoski  et  al.  (1985, 
p.  23). 


A  total  of  nine  space  geodetic  stations  (seven  SLR  stations  and  two  VLBI  stations)  falling 
in  this  regional  vertical  datum  are  considered  in  this  study.  Location  of  these  space  geodetic 
stations  are  shown  in  Figure  4.1.  The  precise  heights  of  these  stations  as  listed  in  the  NASA 
Space  Geodesy  Program  Catalogue  of  Site  Information  (NASA  Tech,  memo  4482, 1993)  refer  to 


44 


the  NGVD29  height  system.  The  VERTCX>N  software  (Version  l.(X),  1992)  received  from 
Vertical  NetwOTk  Branch,  NGS  which  computes  the  modeled  difference  in  the  oithometric  height 
between  the  NAVD88  and  the  NGVD29  for  a  given  location  specified  by  latitude  and  longitude, 
was  used  to  convert  the  NGVD29  normal  orthometric  heights  to  NAVD88  Helmert  orthometric 
heights.  The  amount  of  correction  to  be  applied  at  each  of  the  geodetic  stations  is  shown  in 
Table  4.7.  The  3'x3'  gravity  anomaly  ^d  for  the  conterminous  US,  used  as  input  data  set  for  the 
development  of  GEOID90  Geoid  Height  Model  by  NGS  was  also  used  in  this  study.  These 
gravity  values  are  atmospherically  corrected  Helntert  anomalies  based  on  the  International 
Gravity  Standardization  Net  1971  (IGSN71).  Detailed  documentation  on  this  gravity  data  set 
can  be  found  in  Milbert  (1991a).  As  quoted  by  Milbert  (1991b)  these  3'x3'  mean  anomalies  have 
a  standard  deviation  of  ±1 .5  mgal. 

4.1.2  Ordnance  Datum  Newlvn.  ODN  (England! 

ODN  has  been  established  based  on  the  ddal  observations  taken  between  1st  May,  1915 
to  30th  April  1921  at  the  tide  gauge  on  the  South  pier  of  Newlyn,  Cornwall,  England  (Wilson, 
J.I.,  personal  communication,  March  1993).  Heights  are  referred  to  the  Tide  Gauge  Bench  Mark 
which  is  defined  to  be  4.751  meters  above  the  datum.  Heights  determined  from  this  datum  are 
based  on  either  the  2nd  Geodetic  leveling  of  England  and  Wales  (1920  to  1952)  or  the  3rd 
Geodetic  leveling  of  England  and  Wales  (1953  to  the  present). 

The  SLR  station  (CDP  No.  7840)  at  Herstmonceux,  England  is  the  lone  station  falling  in 
this  regional  datum  to  be  used  in  this  study.  The  station's  oithometric  height  referring  to  the  SLR 
axis  has  been  accurately  determined  in  1991  by  first  order  leveling. 

The  gravity  data  used  around  this  station  are  taken  from  the  same  6'xlO'  atmospherically 
corrected  gravity  data  set  used  by  Despotakis  (1987)  in  his  dissertation  work,  and  they  are  not 
corrected  for  terrain. 

4.1.3  Institut  G^ographique  National  69.  IGN69  (France) 

The  SLR  station  (CDP  No.  7835)  of  Grasse,  France  is  the  only  station  in  IGN69  used  in 
this  study.  IGN69  is  a  normal  height  system  using  <3RS67,  with  zero  at  the  Marseille  tide  gauge. 
The  height  of  the  axis  of  the  SLR  station  used  resulted  from  an  adjustment  of  the  1st  order 
leveling  net  re-observed  during  1961-1969  (C.  Boucher,  personal  communication,  April  1993). 

The  6'xlO'  atmospherically  corrected  (but  not  corrected  for  terrain)  data  set  used  by 
Despotakis  (1987)  was  also  used  in  this  study. 

4.1.4  Normal-Null  (NNl.  (Sea  Level  Datum).  Germany 

The  heights  in  the  German  height  system  (NN)  are  referred  to  the  level  surface  through 
the  reference  point  of  mean  sea  level,  lliis  reference  point  is  situated  37.000m  below  the 
standard  benchmark  which  was  established  via  leveling  lines  from  the  tide  gauge  at  Amsterdam 
("Normaal  Amsterdamsch  Peil,  NAP).  Details  regarding  this  height  system  are  discussed  in 
Torge  (1991,  Section  6.2.3,  p.  228). 

The  SLR  station  (CDP  No.  7834)  at  Wettzell,  Germany  is  used  in  this  study  for  realizing 
the  GVD.  The  gravity  data  used  are  from  the  6'xIO'  mean  anomaly  data  set  used  by  Despotakis 
(1987). 


45 


K  iTTi*  •  wR  IT*  •JiroTiIlf J  ■  fil  •  •  •j J  IWiTiTTiTi 


The  Australian  Height  Datum  was  established  on  May  5. 1971  by  the  division  of  National 
Mapping,  Australia  (National  Mapping  Council,  1979)  canying  out  a  simultaneous  adjustment  of 
97320  km  two-way  leveling,  holttog  mean  sea  level  for  19^^  fixed  at  zero  at  thirty  tide  gauge 
stations  around  the  coast  of  the  Australian  continent.  The  important  characteristics  of  this  &tum 
as  discussed  by  Morgan  (1992)  and  Rizos  et  al.  (1991)  are  listed  as  follows: 

•  ADH71  is  primarily  a  hom^eneous  third  order  or  mapping  datum  with  an  estimated  internal 
precision  of  about  ±8  mm-y^  where  k  is  the  length  of  level  run  in  kilometers. 

•  Data  making  up  the  AHD71  are  not  accurate  enough  for  geodetic,  geodynamic  and 
oceanographic  studies.  A  mid  seventies  re-leveling  program  indicates  that  substantial  errors 
exist  in  the  original  data. 

•  The  difference  between  the  constrained  solution  AHD71  and  a  free  adjustment  suggests  a  rise 
of  mean  sea  level  up  the  east  coast  of  about  l.Sm,  though  this  result  is  in  disagreement  with 
physical  oceanographic  measurements  (Sturges,  1974). 

Two  SLR  stations,  one  at  Yairagadee  (CDP  No.  7090)  in  Western  Australia  and  another 
at  OiToral  VtQley,  Canberra  ((TDP  No.  7943)  in  Eastern  Austrdia  are  considered  in  this  study  for 
realizing  the  Figure  4.3,  shows  the  locations  of  these  stations.  The  gravity  data  used  are 
from  the  2'x2'  mean  anomaly  data  used  by  Despotakis  (1987).  Again  these  data  also  are 
corrected  only  for  atmosphere  and  not  for  teirain. 

4.1.6  Scandinavian  Datum 

Though  there  is  no  single  regional  vertical  datum  known  as  the  Scandinavian  datum,  the 
individual  datums  RH70  (Sweden),  N60  (Finland),  and  NN1954,  NNN1957  (Norway)  were 
combined  to  serve  as  a  single  regional  datum  in  this  study  mainly  for  computational 
convenience.  The  following  factors  also  prompted  the  decision: 

•  These  three  datums  are  in  the  same  peninsular  region  and  also  have  been  included  in  the 
UELN-73  adjustment. 

•  In  the  three  systems,  the  elevations  can,  in  principle  be  reduced  to  any  epoch  by  using 
published  land  uplift  rates  (Ooss  et  al.,  1987,  p.  207). 

•  Most  of  the  ongoing  projects  of  scientific  interest  in  the  Nordic  area  consider  these  datums 
together.  Fennoscandian  Land  uplift  studies  (Ekman  and  MSkinen,  1990),  the  Baltic  sea 
level  project  with  GPS  (Pan  and  Sjdberg,  1993)  are  some  examples  for  the  same. 

•  We  have  a  precise  regional  geoid  model  NKG-89  (Forsberg,  1990)  covering  all  these  three 
datums. 

The  VLBI  site  (CDP  No.  7601)  of  Metsahovi,  Finland,  the  VLBI  station  (CDP  No.  7602) 
at  Tromsd,  Norway  and  the  GPS  station  identified  as  NDG  (Sweden)  in  the  European  north- 
south  GPS  traverse  proposed  by  lAG  SSG  3.88  to  control  and  improve  the  European  geoid 
(Torge  et  al.,  1989)  are  the  three  space  geodetic  stations  considered  in  this  study.  The  ellipsoidal 
coordinates  referring  to  WGS84  and  the  normal  height  referring  to  UELN-73  for  the  GPS  station 
in  Sweden  were  obtained  from  the  Appendix  to  the  paper  by  Torge  et  al.  (1989).  The  precise 
heights  of  the  VLBI  stations  at  Metsahovi  and  Troms^  referring  to  their  local  height  systems 
were  obtained  from  NASA  Tech  Memo  4482.  Since  the  GPS  station  height  is  referr^  to  UELN- 
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73.  and  in  general  Swedish  hei|^ts  are  qjpioximately  equal  to  UELN  heights  (widiin  1-2  cm) 
(Forsberg,  personal  conununication.  May  1993)  the  odier  station  heights  were  also  convened  to 
Swedish  heights  using  the  following  conversion: 

•  Troms^  Norway  local  height  convened  to  Swedish  height  by  subtracting  9  cm  (Forsberg, 
personal  communication,  May  1993). 

•  Metsahovi,  Finland  local  height  to  Swedish  height  by  subtracting  8  cm  (Pan  and  Sjdberg, 
1993). 

The  Faye  anomaly  data  in  a  5^x3°  area  around  each  station  in  Norway,  Sweeten  and 
Finland  were  received  from  Rend  Forsberg.  They  were  point  anomaly  data  sets  selected  closest 
to  the  center  of  every  3'x6'  cell  in  the  area.  Using  the  GECXjRID  program,  in  collocation  mode, 
and  after  removing  the  reference  field  using  OSU91 A  geopotential  model,  the  anomaly  data  were 
predicted  at  every  2*x2'  grid  corner.  The  reference  anomaly  at  the  grid  comers  was  then  added  to 
the  residual  anomaly  pr^cted  and  the  final  data  set  was  used  in  the  computations. 
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MAP  OF  GEODETIC  SITES  -  UNITED  STATES 
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Location  of  Space  Geodetic  Stations  in  U.S.A.  Used  in  This  Study  (NAVD88  Datum) 


LOCATION  OF  SPACE  GEODETIC  STATIONS  -  EUROPE 


2  Location  of  Space  Geodetic  Stations  in  Europe  (covering  Scandinavian.  ODN, 
1GN69,  NN  regional  vertical  datums)  Used  in  This  Study 
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STATIONS  USED  IN  VERTICAL  DATUM  DEFINITION 


soruiiu') 


Figure  4.4  Location  of  Space  Geodetic  Stations  Used  in  Global  Vertical  Datum  Definition 


Table  4. 1  StadcMis  Used  in  Global  Vertical  Datum  Definition  Study 


Station  Latitudo  Longitude  Ellip.ht  Elevation 

(a-6378I36.3)  above 

Local  OatuB 
(■)  (B) 


7051  QUINCY 

39 

58 

24.577 

239 

3 

37.432 

1060.686 

1083.480 

7062  SANDOI 

32 

36 

2.667 

243 

9 

32.660 

989.219 

1022.000 

7082  BEARLX 

41 

56 

0.899 

248 

34 

45.417 

1963.700 

1976.515 

7091  WESTPD 

42 

37 

21 i 701 

288 

30 

44.233 

92.794 

120.502 

7086  FTDAVS 

30 

40 

37.311 

255 

59 

2.369 

1962.178 

1983.161 

7105  GSFCMD 

39 

1 

14.176 

283 

10 

20.050 

19.890 

53.050 

7204  GREENB 

38 

26 

16.164 

280 

9 

50.944 

813.222 

844.715 

7110  MOPSAK 

32 

53 

30.254 

243 

34 

38.140 

1839.702 

1870.778 

7234  PITOWN 

34 

18 

3.666 

251 

52 

50.665 

2365.390 

2385.940 

7069  PATRIK 

28 

13 

40.669 

279 

23 

39.191 

-22.937 

6.533 

7602  TROMSO 

69 

39 

46.702 

18 

56 

21.964 

133.599 

101.400 

7601  NETSAH 

60 

14 

31.069 

24 

23 

2.999 

60.269 

40.598 

1001  ONSGPS 

57 

53 

50.713 

12 

2 

46.382 

39.959 

2.736 

7835  GRASSE 

43 

45 

16.883 

6 

55 

16.034 

1323.560 

1272.100 

7834  WETZLL 

49 

8 

41.775 

12 

52 

41.142 

661.836 

614.440 

7840  RGOUKl 

50 

52 

2.565 

0 

20 

10.040 

76.098 

30.693 

7943  CANBRR 

-35 

-37 

-29.751 

148 

57 

17.296 

949.635 

929.500 

7090  YARRAG 

-29 

-2 

-47.410 

115 

20 

48.272 

242.051 

266.568 
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4.2  Systematic  Coircctions  ip  Gravity  Anomalies 

The  systematic  conections  needed  to  be  applied  to  the  terrestrial  ^vity  anomalies  were 
mainly  the  horizontal  datum  inconsistency  correction  (5gH)  and  the  ^vity  fmmula  correction 
(Sgc).  as  they  have  been  already  corrected  for  atmospheric  effects.  Since  we  didn't  have  dense 
DTM  data  for  applying  the  terrain  correction  to  the  gravity  anonutiy  data  sets  that  were  not 
corrected  for  the  terrain,  they  were  used  as  such.  Despotakis  (1989),  based  on  his  undulation 
computations  with  data  at  12  space  stations  in  the  western  United  States,  has  shown  that  the 
mean  and  RMS  difference  of  computed  undulations  decreased  by  about  SO  cm  when  the  terrain 
corrections  were  taken  into  account.  The  standard  deviation  decreased  about  20  cm  when 
terrain-corrected  gravity  data  was  used.  The  extent  of  5gH  and  Sge  corrections  needed/applied  to 
the  various  gravity  data  sets  are  discussed  below: 

4.2.1  Horizontal  Datum  inconsistanev  Collection  (Sch) 

As  mentioned  in  section  3.1.1(b),  this  correction  is  needed  because  of  the  systematic  error 
introduced  in  the  computation  of  normal  gravity  due  to  erroneous  latitude  values.  Detailed 
theory  behind  such  a  correction  is  given  in  Heck  (1990,  Section  3.3,  pp.  96-101).  In  Table  4.2, 
the  (tetails  regarding  the  size  and  shape  of  the  reference  ellipsoids  associated  with  the  local 
geodetic  datums  are  given.  Also  changes  in  the  equatorial  r^us  Aa  and  in  the  flattening  Af 
between  the  locally  adopted  reference  ellipsoids  and  WGS84  reference  ellipsoid  are  given.  Table 

4.3  gives  the  datum  transformation  parameters  from  the  various  local  horizontal  datums  on  which 
the  stations  used  in  this  study  lie  and  the  ITRF91  reference  frame.  The  transformation  is  done  in 
two  stages.  In  the  first  stage,  the  transformation  is  done  from  Local  geodetic  datum  to  WGS84 
datum  using  the  parameters  given  in  Table  7.13  in  DMA  Technical  Report  Part  1  and  in  the 
second  stage  WGS84  to  mu^l  combining  the  parameters  given  in  l^ucher  and  Altamimi 
(1991,  Table  4,  p.  53)  and  lERS  Tech.  Note  12  (Boucher  et  al.,  1992a,  Table  2,  p.  27). 

Table  4.2  Local  Geodetic  System  -  To  -  WGS84;  Comparison  of  Ellipsoidal  Parameters 
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Table  4.3  Hivizontal  Datum  Transformation  Parameters  from  Local  Geodetic  Datum  to  ITRF91 


The  effect  of  horizontal  datum  inconsistencies  on  the  computed  gravimetric  quantities  are 
given  in  Table  4.4.  Since  the  latitude  values  of  the  stations  in  the  same  regional  vertical  datums 
v^  considerably,  the  extent  of  correction  also  changes  from  one  station  to  another.  The  details 
given  in  Table  4.4  were  computed  only  fra-  one  station  in  each  of  the  different  datums  considered. 
The  name  of  the  station  at  which  the  computation  was  carried  out  is  given  in  parenthesis. 

Table  4.4  Errors  Introduced  due  to  Horizontal  Datum  Inconsistencies  on  Certain  Gravimetric 
Quantities 


From  the  table  above,  it  is  clear  that  the  correction  required  due  to  horizontal  datum 
inconsistencies  is  quite  substantial  and  must  be  taken  into  account  when  ±10  cm  accuracy  is 
sought  in  the  definition  of  the  GVD.  The  computed  values  in  Table  4.4  compare  well  with  the 
values  given  in  Tables  2a  and  4  of  Heck  (1990).  The  gravity  data  used  in  this  study  were  all 
conect^  for  the  Horizontal  datum  inconsistencies  (Sgu). 


The  gravity  formula  used  in  the  computation  of  gravity  anmnalies  in  NAVD88  and  in  the 
Scandinavian  datums  refer  to  the  formula  implied  by  the  GRS80  constants,  in  which  the  semi¬ 
major  axis  is  taken  as  6378137.0m.  Similarly  in  other  European  and  Australian  stations,  as 
discussed  in  Despotakis  (1987),  the  gravity  formula  is  computed  using  a  =  6378136.(hn.  But  in 
our  study  for  the  adopted  ideal  reference  ellipsoid  we  have  a  =  6378136.3m.  Therefore,  the 
gravity  anomalies  are  to  be  corrected  for  change  in  semi-major  axis  using  the  procedure 
explained  in  Despotakis  (ibid.,  p.  51).  The  computi^  5gG  corrections  for  the  gravity  data  sets  in 
NAVD88  and  the  Scandinavian  datums  amounts  to  -0.216  mgal  and  for  other  gravity  data  sets  in 
ODN,  IGN69,  NN  and  AHD71  datums,  it  is  +0.093  mgal. 
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Two  different  techniques,  namely.  Modified  Stokes'  technique  and  Least  Squares 
Collocation  (LSC)  technique  were  used  in  the  computation  of  gravimetric  height 
anomaly^ndulation  at  the  space  geodetic  stations  used  in  the  GVD  definition  study. 

4.3. 1  Modified  Stokes'  Technique 

As  explained  in  section  3.3.1,  this  technique  makes  use  of  tenestrial  gravity  anonudy  data 
in  a  cap  around  the  station  and  a  set  of  potential  coefficients  from  an  accurate  geopotential  model 
to  compute  the  height  anomaly/undulation.  Equations  (3.71)  and  (3.72)  are  the  key  equations 
used  for  this  purpose. 

4.3. 1 . 1  Selection  of  Cap  Size 

Accuracy  of  computed  height  anomaly/undulation  depends  largely  on  the  extent  of  the 
tenestrial  gravity  data  used  around  the  computation  point  As  shown  in  Table  8  of  Pavlis  (1991, 
p.  1 19),  the  grid  spacing  has  only  a  minor  effiect  on  the  quali^  of  the  result  Computation  of 
undulation^eight  anomaly  using  Modified  Stokes'  technique  with  (3.71)  or  (3.72),  involves  four 
main  error  sources  in  its  estimation.  They  are: 

1.  Errors  associated  with  the  gravity  anomalies  being  given  at  discrete  locations  (or  as  mean 
values)  instead  of  a  continuous  function  (called  sampling  error  or  discretion  error  and 
denoted  as  5Ni). 

2.  Errors  associated  with  gravity  data  in  the  cap  surrounding  die  station  (called  propagated  error 
and  denoted  as  SN2). 

3.  Errors  associated  with  the  potential  coefficients  of  the  reference  field  (called  commission 
error  and  denoted  as  5N3). 

4.  Errors  associated  with  the  neglected  potential  coefficients  above  degree  M  (called  omission 
error  and  denoted  as  5N4). 

The  total  global  root  mean  square  error  assuming  uncorrelated  error  sources  is: 

5N  =  (SnJ  +  8n5  +  5N§  +  8n3)1/2  (4.2) 

Table  4.5  below  shows  the  estimated  accuracy  of  gravimetric  undulation  computed  using 
equation  (4.2)  for  different  cap  sizes  of  terrestrial  gravity  anomaly  data.  Detailed  discussions  on 
the  different  error  sources  listed  above  and  the  required  expressions  for  computing  the  same  are 
given  in  Section  3.1  of  Pavlis  (1991)  and  in  Section  2  of  Despotakis  (1987).  With  potential 
coefficients  frtnn  the  OSU91 A  model  and  anomaly  errcv  degree  variances  referring  to  mean  earth 
radius  of  6371  kilotneters,  the  values  in  Table  4.S  were  computed  using  the  reciprocal  distance 
error  covariance  nnodel,  described  in  Pavlis  (ibid..  Section  3.4.2),  with  variance  4  mgal^  and 
correlation  length  0?1.  As  can  be  seen  from  the  results,  due  to  the  increase  in  the  propagated 
error  with  an  increase  in  cap  size,  beyond  certain  cap  radius  use  of  terrestrial  gravity  data  does 
not  improve  but  degrades  the  result.  'The  estimated  total  error  in  the  computed  undulation  for  cap 
size  \|/  =  2**  and  y  =  3**  are  almost  the  same  indicating  the  use  of  Ag  data  in  a  T*  cap  size  for 
computing  the  undulations  should  be  adequate. 
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Table  4.5  Estimated  Accuracy  of  Computed  Gravimetric  Undulations  Using  Teirestrial  Gravity 
Data  in  Different  Cap  Sizes  Along  with  the  OSU91A  Global  Geopotential  Model 
(degree  360).  Units  are  in  cm. 


As  mentioned  in  the  previous  section,  the  gravity  data  around  the  space  geodetic  stations 
in  the  Scandinavian  datum  were  available  only  in  an  area  of  in  the  latitude  and  longitudinal 
directions.  Since  these  stations  are  located  in  extreme  northern  latitudes,  the  gravity  data 
required  in  the  longitudinal  direction  is  almost  double  the  latitude  extent,  which  can  be  computed 
using  the  following  relation: 

AX  =  -^  (4.3) 

cos<|> 

where  0  is  the  latitude  of  the  station,  A0  and  AX  are  the  latitude  and  longitudinal  extents.  As  the 
data  availability  in  longitudinal  direction  from  the  computation  point  is  limited  to  2**  to  2?5, 
terrestrial  gravity  data  to  a  smaller  latitude  extent  alone  can  be  used  Table  4.6  shows  the  extent 
of  the  cap  size  radius  of  the  terrestrial  gravity  data  that  can  be  used  for  estimating  height 
anomaly/undulation  with  Modified  Stokes'  technique,  with  available  data. 

Table  4.6  CTomputation  of  Possible  (Hap  Sizes  of  Terrestrial  Gravity  Data  That  (Han  Be  Used  in 
Scandinavian  Datum 


tations 


Latitude  0 


Onsala  (GPS), 
Sweden 


Tables  4.7  and  4.9  show  the  results  of  the  ^vimetric  height  anomaly/undulation 
computations  carried  out  at  the  space  geodetic  stations  in  the  United  States,  Europe  and 
Australia.  Since  ODN  (England)  values  are  given  as  orthometric  heights,  the  gravimetric 
undulation  computed  at  the  Herstmonceux  SLR  station  (CDP  No.  7840)  alone  is  shown 
separately  in  Table  4.9(b).  For  other  European  stations  the  results  are  given  in  Table  4.9(a). 

The  first  column  in  all  the  tables  refers  to  the  number  and  name  of  the  stations  and  the 
next  five  columns  show  the  various  components  of  the  height  anomaly/undulation  computed 
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using  equaticms  (3.71)  and  (3.72).  The  various  components  shown  in  column  (2)  to  (6)  are:  the 
cap  contribution  in  (2),  remote  zone  contribution  computed  using  potential  coefficients  in  (3), 
indirect  effect  or  corrections  to  height  anomaly  dejtending  on  computed  quantity  Nj  otCI  in  (4). 
ellipsoidal  corrections  in  (5)  and  the  total  gravimetric  undulation/height  anomaly  in  (6). 

In  Table  4.7,  the  geometric  undulation  at  the  U.S.  stations  computed  with  NGVD29 
heights,  corrections  to  conven  these  undulations  to  NAVD88  height  system  computed  using 
VERT(jON  (version  1.00)  software  and  the  resulting  undulations  in  NAVD88  height  systems  are 
also  given  in  columns  (7),  (8)  and  (9)  respectively.  Similarly,  in  Table  4.9(a),  the  geometric 
height  anomaly  comput^  with  respect  to  the  local  height  datum,  corrections  needed  to  convert  to 
Swedish  height  systems  in  the  case  of  station  numbers  1001  and  7602  and  the  final  height 
anomaly  are  given  in  columns  (7)  to  (9).  Since  the  other  European  stations  and  Australian 
stations  are  assumed  to  lie  in  their  local  height  datum,  no  correction  was  needed  to  the  computed 
geometric  height  anomaly/undulation  at  these  stations. 
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Since  the  LSC  technique  is  based  on  the  ergodic  assumption  discussed  in  Section  3.3.2, 
the  basic  covariance  functions  referred  to  in  equations  (3.75)  and  (3.76)  are  assumed  to  be 
homogeneous  and  isotropic,  depending  only  on  the  spherical  distance  between  two  points.  In 
this  study  theoretical  covariances  were  computed  baskl  on  the  noise  model  of  the  reference  field 
and  on  a  hypothetical  model  for  the  behavior  of  the  anomaly  degree  variances  (Cn)  of  the  Earth's 
gravity  field.  The  expressions  for  the  fundanrental  covariance  functions  used  in  this  study,  as 
given  in  Ba]$i6  and  Rapp  (1992),  are  as  follows: 


360. 


C(Ag,Ag)=  l8C„s“*‘P„(cosv)+  IC„s“"*P.(cos\|r) 

n-2  361 


C(Ag,N)=-^ 

Y 


360  ^ 

2:^s"*‘P„(cos  v)+ 2:-bi^s”*"P„(cos  V) 
ii-2n  - 1  36!  n  - 1 


C(N,N)  =  i^ 


[n-2(n-iy 


•s“"‘P„(cosv)+I 


5;(n-l)" 


s"^  P„(cos\p) 


(4.4) 

(4.5) 

(4.6) 


where  SCn  are  the  error  degree  variances  related  to  the  reference  field  defined  through  a  set  of 
geopotentid  coefficients  computed  as. 


5C„  = 


GM 


in-lf 


\2(n+2) 


S(e; 


4.?* 
cam  '  *'mr 


ibkO 


(4.7) 


^cnmt'^snin  ^  Standard  deviations  of  the  potential  coefficients  Cnm  and  S  nm>  Re  is 
the  mean  es^  ra^us  (6371  km),  Rb  is  the  radius  of  Bjerhamar's  sphere  and  n,  m  are  degree  and 
order  of  the  harmonics, 

Cn  are  the  anomaly  degree  variances  computed  by  the  well  known  Tscheming/Rapp 
(1974)  (T/R)  anomaly  degree  variance  model. 


C  -  A(n-l) 

"  (n-2)(n  +  B) 


(4.8) 


with  parameters  A  =  343.3408  mgal2,  B  =  24  and  s  =  J  =  0-9988961  recommended  by 

Jekeli  (1978). 

The  covariance  functions  computed  with  the  error  degree  variances  of  the  OSU91A 
reference  field  up  to  degree  360  and  beyond  in  the  T/R  model  with  Jekeli  parameters  and  used  by 
Ba]^i6  and  Rapp  (1992)  were  also  used  in  this  study.  The  behavior  of  the  covariance  functions  are 
shown  in  (ibid..  Table  4.1  and  Figure  4.4).  For  the  prediction  of  residual  geoid  undulation/height 
anomaly  using  the  local  gravity  data  around  the  computation  point,  the  global  covariance 
function  should  be  scaled  to  the  local  covariance  function  so  as  to  reflect  appropriately  the 
accuracies  of  predicted  quantities. 

The  scaling  procedure  described  below  follows  closely  the  discussion  in  Hwang  (1989, 
Section  3.3). 
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For  the  predictim  area,  the  variance  of  residual  anomalies  is  conqtuted  by. 


(4.9) 


where  n  is  the  total  number  of  points  used  for  prediction  and  Ag  is  the  mean  of  the  residual 
anomalies. 


Let  CXjGCO)  be  the  variance  of  anomalies  based  on  (4.4).  We  have  to  tq>ply  a  scaling  factor  a  to 
this  quantity  so  that  the  resultant  will  be  identical  to  the  variance  of  anomalies  given  by  equation 
(4.9).  Thus  a  is: 


a  = 


CX3G(0) 


(4.10) 


The  local  covariance  functions  then  can  be  constructed  by  multiplying  equations  (4.4),  (4.5)  and 
(4.6)  by  the  scaling  factor  a . 


To  ensure  a  realistic  scaling  factor,  the  lowest  value  of  a  was  taken  to  be  0.120  and  the 
maximum  value  limited  to  1.5.  In  this  range,  it  was  assumed  that  the  accuracy  of  predicted 
undulation  will  be  within  acceptable  limits. 

Based  on  the  above  discussions,  the  equations  (3.75)  and  (3.76)  are  modified  for  this 
study  as  follows: 


S=C.[c.+ic_]''Y 
E.=o|c.-C.[c.+ic„]  C^] 


(4.11) 

(4.12) 


Equations  (4.11)  and  (4.12)  are  similar  to  the  expressions  given  by  Hwang  (1989)  as  equations 
(3.24)  and  (3.25). 

4.3.2.2  Numerical  Results 

The  residual  undulations  computed  at  the  space  geodetic  stations  in  the  United  States  are 
shown  in  Table  4. 10.  The  conqiutations  were  carried  out  with  different  grid  sizes  and  number  of 
points  to  understand  their  effects  on  computed  undulations.  Residual  undulations  computed  with 
gravity  data  in  6'x6'  grid  size  and  covering  an  area  of  4'*x4‘’  around  the  computation  point  were 
select^  for  this  study.  This  was  done  not  only  fmr  the  reason  that  it  covered  the  same  exrent  of 
gravity  data  used  by  the  Modified  Stokes'  technique  to  compute  the  height  anomaly/undulation, 
but  of  the  three  cases  considered  they  had  the  least  standard  deviations  in  their  estimates.  Table 
4.11  shows  the  estimated  height  anomaly  along  with  the  accuracy  estimate  at  the  Scandinavian 
datum  stations  as  well  as  at  the  other  European  space  geodetic  stations  considered  in  the  study. 
The  number  of  points  used  in  predicting  the  residual  height  anomaly/undulation  are  shown 
within  the  parenthesis  under  each  computed  value.  The  undulation  computed  at  the  Australian 
stations  are  shown  in  Table  4.12. 
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Table  4.10  G>mputation  of  Residual  Undulation  with  Different  Grid  Sizes  and  Total  Undulation  at  the  U.S.  Space  Geodetic  Stations 


Table  4. 1 1  Computation  of  Height  Anomaly/Undulation  at  Scandinavian  and  Other  European 
Stations  Using  LSC  Technique.  Units  are  in  meters. 


Table  4.12  Computation  of  Gravimetric  Undulation  at  Ausoidian  Space  Geodetic  Stations  Using 
LSC  Technique.  Units  are  in  meters. 


Tables  4.13  to  4.15  show  the  comparison  between  the  gravimetric  undulation/  height 
anomaly  computed  using  the  Modified  Srokes'  technique  and  the  LSC  technique  with  the 
geometric  undulations  computed  at  the  stations  for  the  idi^  reference  ellipsoid  (a  =  6378136.3 
m).  The  undulations  interpolated  from  regional  geoid  height  models  GEOID93  (in  die  United 
States),  NKG-89  (in  Scanchnavian  countries)  and  AUSGEOID93  in  Australia  received  from  Jim 
Steed  (personal  communication,  Dec.  93)  are  also  included  in  the  tables  for  better  comparison. 
The  values  of  the  observation  vector  Y  in  equations  (3.41)  and  (3.63)  are  also  shown  for  the 
different  techniques  with  corresponding  statistical  information  such  as 
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iDean  =  x 


r.m.s.  values 


s  tan  danl  deviation  = 


(4.13) 
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The  results  from  Tables  4.13  and  4.14(a)  and  (b)  show  that  the  computation  of 
gravimetric  undulation/height  anomaly  using  the  two  different  techniques,  i.e..  Modified  Stokes' 
and  LSC  techniques,  are  consistent  at  ±20  cm,  except  for  the  four  stadms  7082, 7204, 7086  and 
7069  falling  in  Ae  United  States,  where  the  difference  between  the  two  techniques  varies  from 
29  cm  to  54  cm.  In  Table  4.15  also,  we  see  the  undulation  computed  at  the  two  Australian 
stations  using  the  Modified  Stokes'  and  LSC  techniques  differ  by  about  45  cm  to  65  cm.  One  of 
the  reasons  for  such  a  large  difference  between  the  computed  gravimetric  undulations  from  the 
two  techniques  may  be  that  only  a  limited  number  of  data  points  (1600)  was  used  in  the  LS(r 
technique  for  predicting  the  residual  undulation  due  to  computer  limitaticxis  in  memory  and  CPU 
time.  The  following  Table  4.16  show  the  gravimetric  undulations  computed  at  the  SLR  station 
7943  Orroral  Valley,  Canberra,  Australia  using  Modified  Stokes'  technique  with  terrestrial 
gravity  data  in  different  cap  sizes,  the  undulation  value  computed  using  the  LSC  technique  with 
the  data  points  limited  to  1600  and  also  the  geoid  undulation  value  interpolated  from 
AUSGEOID93  model  received  from  Jim  Steed  (personal  communication,  Dec.  1993). 

Table  4.16  Computation  of  Gravimetric  Undulation  Using  Modified  Stokes'  Technique  With 
Terrestrial  Gravity  Data  in  Different  (}ap  Sizes.  (Units  are  in  meters) 


("‘with  data  points  for  prediction  limited  to  1600  at  Orroral  Valley,  Canberra) 


From  Table  4.16,  we  see  that  the  undulations  computed  using  the  terrestrial  gravity  data  in  a  cap 
size  of  0?5,  can  vary  to  an  extent  of  about  40  cm,  when  compared  with  undulations  computed 
with  gravity  data  in  a  cap  size  of  2?0.  The  geoidal  undulations  obtained  from  the  AUSGEOID93 
model  which  uses  terrestrial  gravity  data  in  a  cap  size  of  0^5  and  OSU91A  potential  coefficients 
for  its  computation  (Kearsley  et  al.,  1990)  compares  well  with  the  undulations  computed  using 
the  LSC  technique  with  the  number  of  data  points  used  for  prediction  limited  to  1600  and  also 
with  the  undulation  computed  using  the  Mo^fied  Stokes'  technique  with  cap  size  =  0?5.  This 
suggests  that  the  differences  in  the  computation  of  undulation  between  the  Modified  Stokes'  and 
LSC  techniques  may  be  due  to  the  extent  of  the  terrestrial  gravity  data  used  in  the 
computation/prediction. 

In  Table  4.15  we  also  notice  a  difference  of  about  1.0  meter  between  the  computed  Y 
values  at  the  Yarragedee  station  in  Western  Australia  and  Orroral  Valley  station  in  eastern 
Australia.  The  difference  observed  is  for  both  the  Modified  Stokes'  and  LSC  techniques  used  for 
computing  the  gravimetric  undulations,  as  well  as  for  the  AUSGEOID93  undulation  values 
received  from  Jim  Steed  (personal  communication,  Dec.  93),  though  both  the  stations  are 
assumed  to  lie  in  the  same  height  system  AHD7i.  Since  as  quoted  by  J.  Steed  (ibid.), 
AUSGEOID93  is  the  most  accurate  geoid  model  available  in  Australia  at  present  and  also 
ellipsoidal  heights  of  the  stations  are  obtained  to  an  accuracy  of  ±5  cm,  the  large  difference  of 
1.0  meter  in  computed  Y  value  (Y  =  h  -  HahD71  -  Ng)  between  the  two  stations  as  observed  from 
all  the  three  techniques,  causes  serious  concern.  This  can  be  due  to  one  or  more  of  the  following 
factors: 
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(a)  the  long  wavelength  error  in  the  ^obal  geof^ntial  model  used, 

(b)  bias  caused  by  systematic  errors  in  ^  gravity  anomaly  data  used,  (for  example  the 
gpvity  anomaly  around  the  Australian  stations  are  not  corrected  ftn*  terrain),  and 

(c)  distortions  in  the  defined  height  system. 

To  verify  whether  the  difference  in  Y- value  is  mainly  due  to  long  wave  length  errors  in 
the  OSU91A  global  geopotential  model  used  in  the  study,  the  computations  of  the  height 
anomaly/undulation  using  ModiHed  Stokes'  technique  were  repeated  with  the  latest  JGM'2 
(NASA  model)  from  degree  2  to  70  and  augmented  with  OSU91A  model  from  degree  71  to  360 
and  the  results  are  shown  in  Appendix  'A*  of  this  report.  Details  regarding  JGM-2  model 
development  can  be  found  in  Nerem  et  al.  (1993).  From  the  results,  we  see  the  difference  in  the 
value  of  Y  computed  using  the  different  potential  models  at  the  Australian  space  geodetic 
stations  vary  only  by  5  to  9  cm,  implying  that  die  long  wavelength  error  in  die  geopotential 
model  used  alone  is  not  the  contributor  to  the  1.0  m  difference.  Also  the  terrain  corrections 
computed  at  the  Australian  stations  using  ETOPOSU  S'xS'  elevation  data  set  (as  no  other  dense 
elevation  data  are  available  in  this  area),  shows  the  terrain  correction  is  on  the  order  of  0.07 
mgals  accounting  for  an  undulation  correction  of  about  2  cm.  The  above  discussions  point  out 
that  the  large  difference  in  computed  Y- values  between  the  two  stations  may  be  mainly  due  to 
the  distortions  in  the  AHD71  height  system.  As  discussed  by  Morgan  (1992)  and  referred  to  in 
Section  4.1.5,  there  is  reason  to  believe  that  substantial  eirors  must  exist  in  the  original  leveling 
data  used  for  developing  the  AHD71  height  system. 


The  estimated  values  of  parameters  defining  the  GVD,  that  is  AWq  and  C^QjO  (ft>f  »  *  1. 
2,...n,  n  being  the  number  of  regional  vertical  datums  considered)  and  their  accuracy  estimates 
were  computed  using  equations  (3.87)  and  (3.88)  after  setting  up  the  observation  equations  as 
explained  in  Chapter  in  at  each  of  the  17  stations  considered  in  this  study.  The  elements  of  the 
A  matrix  (of  size  17x7)  used  in  these  equations  are  shown  in  Table  4.17,  where  its  first  column 
correspond  to  the  AWq  term  and  the  rest  six  columns  to  the  C^o  values.  The  elements  of  the  K 
vector  were  chosen  to  be  equal  to  the  number  of  space  geodetic  stations  used  in  the 
corresponding  regional  vertical  datum.  As  stated  in  Section  3.4,  the  value  in  K  corresponding  to 
the  AWo  parameter  was  chosen  to  be  equal  to  zero. 


The  variance-covariance  matrix  of  observations  £y  ^as  set  up  assuming  that  the  three 
observations  types:  ellipsoidal  height  (h),  orthometric  height  (H)  and  height  anoi^y/undulation 
(N)  are  independently  determined  such  that; 

2:y  =  Ih  +  +£n  (4.4) 

where  £h>  the  variances  of  the  observed/computed  quantities  h,  H  and  N 

respectively. 

Determination  of  orthometric  height  standard  deviation  is  rather  a  difficult  task. 
Conceptually,  one  would  take  the  diagonal  elements  of  the  variance-covariance  matrix  of  the 
vertic^  network  adjustment  for  a  specific  datum.  Since  such  information  is  not  available  for  all 
the  datums  considered  in  this  study,  based  on  the  general  information  regarding  the  regional 
leveling  networks,  certain  reasonable  standard  deviations  can  be  assumed.  For  the  NAVD88 
datum  heights  obtained  by  applying  corrections  to  NGVD29  heights,  a  standard  deviation  of  ±5 
cm  has  been  assumed  at  all  the  stations  in  the  United  States.  Similarly,  at  the  Australian  stations. 


an  accuracy  of  ±10  cm  has  been  assumed,  since  AHD71  is  based  on  only  a  third  order  leveling 
network  (MOTgan.  1992).  For  the  Eurqsean  stations,  the  standard  deviations  of  the  normal  height 
has  been  interpolated  from  the  flgures  of  the  adjusted  UELN73  leveling  networics  given  in 
Ehmsperger  et  al.  (1982). 

Determination  of  ellipsoidal  height  standard  deviation  at  the  space  geodetic  stations  can 
be  easily  done,  since  Table  T6  in  lERS  Technical  Note  12  (Boucher  et  sQ.,  1992a)  gives  the 
standard  deviations  of  the  rectangular  coordinates  X.  Y,  Z  also,  in  addition  to  providing  station 
coordinates  in  rrRF91  reference  frame.  Assuming  that  the  rectangular  comdinates  are 
uncorrelated,  by  error  propagation  we  can  compute  the  error  in  the  corresponding  ellipsoidal 
coordinates. 

The  determination  of  geoid  undulation  accuracy  can  be  done  based  on  the  discussion  in 
Chapter  5.4  of  Rapp  and  Balasubramania  (1992),  wlKn  using  the  Modified  Stokes'  technique  for 
computing  the  gravimetric  height  anomaly/undulation.  With  OSU91 A  standard  deviations  of  the 
potential  coefficients,  assuming  a  correlation  length  of  0?1  and  using  Gauss-Markov  first  order 
error  covariance  model  (Pavlis,  1991,  p.  107),  the  estimated  stands^  deviations  of  the  height 
anomaly/undulation  are  given  in  Table  4.18.  For  the  LSC  technique  the  standard  deviations  of 
the  pre^cted  undulations  can  be  obtained  from  Tables  4.10  to  4.12.  The  diagonal  elements  of 
the  £y  matrix  for  the  Modified  Stokes'  technique  and  the  LSC  technique  are  shown  in  Tables 
4.19  and  4.20,  respectively. 
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Table  4.17  Elements  of  A  Matrix  With  17  Stations  and  7  Datum  Unknowns  With  Varying  Cap 
Ra^us.  Units  are  in  cm. 


-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

1.0978 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

0.0000 

1.0581 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

0.0000 

1.0390 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

0.0000 

1.0390 

0.0000 

0.0000 

0.0000 

0.0000 

-1.0206 

0.0000 

0.0000 

1.0978 

0.0000 

0.0000 

0.0000 

-1.0206 

0.0000 

0.0000 

0.0000 

1.0978 

0.0000 

0.0000 

-1.0206 

0.0000 

0.0000 

0.0000 

0.0000 

1.0978 

0.0000 

-1.0206 

0.0000 

0.0000 

0.0000 

o.cooo 

0.0000 

1.0978 

-1.0206 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 

1.0978 
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Table  4.19  Computation  of  Variance-Covariance  Matrix  of  Observations  £y  (Modified 
Technique) 


Table  4.20  Computation  of  Variance-Covariance  Matrix  of  Observations  Ly  (Least-Squares 
Collocation  Technique) 
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The  parameters  defining  the  GVD  estimated  separately  for  Modified  Stokes'  technique 
and  LSC  technique  are  shown  in  Figures  4.5  to  4.7.  The  estimated  parameters  are  the  vertical 
separation  between  the  defined  Global  Vertical  Datum  (GVD)  and  the  adopted  reference 
ellipsoid,  AW(/y,  and  the  vertical  separation  between  the  GVD  and  regional  vertical  datum, 
CQfl/y,  for  i  =  1, 2,...n  where  n  is  the  number  of  regional  vertical  datums  considered  in  the  study. 


Figure  4.5  Regional  Vertical  Datum  Separation  From  Defined  Global  Vertical  Datum  (Using 
ModiHed  Stokes'  Technique  to  Compute  N  or  Units  are  in  cm 


Figure  4.6  Regional  Vertical  Datum  separation  from  defined  Global  Vertical  Datum  (using  LSC 
technique  for  computation  of  N  or  0>  units  are  in  cm 


Figure  4.7  Regional  Vertical  Datum  Separation  From  Defined  Global  Vertical  Datum  (From 
LSC  Technique  but  Using  the  Weight  Matrix  of  ModiHed  Stokes'  Technique);  Units 
are  in  cm 


Figures  4.5  and  4.6  show  the  estimated  values  of  the  separation  between  the  regional 
vertical  datums  considered  in  the  study,  as  well  as  the  adopted  reference  ellipsoid  from  tl^.e 
defined  Global  Vertical  Datum  (GVD)  surface  computed  using  Modified  Stokes'  and  LSC 
techniques.  Hiough  from  these  figures  it  is  clear  that  the  different  regional  vertical  datums  used 
in  the  GVD  definition,  show  the  same  pattern  of  their  relative  positions  (except  Scandinavian 


75 


datum),  their  numerical  values  of  separation  are  quite  different.  One  of  the  reasons  for  their 
numerical  difference  is  the  used  weight  matrix,  L~\  for  estimating  the  parameters.  As  can  be 
derived  from  Tables  4.19  and  4.20,  the  estimated  errors  in  the  predicted  height 
anomaly/undulation  using  LSC  technique  is  unbelievably  over-optimistic  at  certain  stations, 
leading  to  a  larger  weight  matrix  than  the  one  used  for  Mo^fied  Stokes'  technique.  Even  though 
the  computed  Y  values  using  Stokes'  and  LSC  techniques  agree  within  ±20  cm  for  most  stations, 
say  particularly  in  Europe,  the  standard  deviations  of  the  computed  height  anomaly/undulation 
using  Stokes'  technique  is  2  to  4  times  larger  than  for  the  predicted  undulation  using  LSC 
technique.  This  is  mainly  because,  the  ergodic  assumption  on  which  the  LSC  technique  is  based, 
does  not  differentiate  between  probabilistic  mean  and  spatial  mean  (which  are  very  different 
when  real  world  ^ta  is  considered)  and  hence  generally  result  in  optimistic  accuracy  estimation. 

Since  the  two  different  techniques  are  mainly  employed  to  compute  the  gravimetric 
height  anomaly/undulation  and  the  results  compare  ve^  well  also  (at  least  at  U.S.  and  European 
stations)  it  should  be  possible  to  use  the  same  variance-covariance  matrix  for  both  the  cases.  The 
main  difference  in  setting  up  the  weight  matrix  for  both  the  cases,  as  can  be  inferred  from  Tables 
4.19  and  4.20,  comes  from  ^e  large  standard  deviations  in  predicted/computed  undulation,  as  we 
have  assumed  the  same  standard  deviations  in  orthometric  heights  and  ellipsoidal  heights  for 
both  techniques. 

Assuming  the  weight  matrix  computed  for  Modified  Stokes'  is  acceptable  to  LSC 
technique  also,  the  parameters  defining  the  GVD  were  estimated  again  using  the  LSC  technique 
with  the  new  weight  matrix.  The  results  are  shown  in  Figure  4.7. 

Figures  4.5  and  4.7  compare  well  except  for  certain  numerical  differences.  Even  the 
Scandinavian  datum  which  came  out  to  be  below  the  German  datum  in  Figure  4.6,  is  placed 
above  in  Figure  4.7  agreeing  with  the  results  in  Figure  4.5.  The  Scandinavian  datum  which  is  58 
cm  above  Ae  defined  GVD  in  Figure  4.5  is  about  39  cm  above  the  same  in  Figure  4.7.  The 
Australian  height  datum  AHD71  which  is  about  102  cm  above  the  GVD  in  Figure  4.5  is  only  52 
cm  in  Figure  4.7.  Tlie  reason  as  can  be  seen  from  Table  4. 15  is  that  the  gravimetric  undulations 
computed  at  the  Australian  stations  using  Modified  Stokes'  technique  and  LSC  technique  differ 
by  about  50  to  70  cm.  The  probable  reasons  for  such  a  difference  are  discussed  in  Section  4.3.4. 
The  datum  differences  between  the  Scandinavian  datum  (assumed  to  coincide  with  Swedish 
Height  System)  and  German  Height  System,  varies  from  4  cm  ±  12  cm  in  Figure  4.7  to  30  cm  ± 
12  cm  in  Figure  4.5.  The  estimates  of  difference  between  the  same  datums  as  given  by  Pan  and 
Sjoberg  (1993)  is  about  16  cm  ±  9  cm.  The  difference  of  33  cm  between  the  IGN69  and  ODN 
height  systems  agrees  well  with  the  estimate  of  30  cm  ±  8  cm  given  by  Willis  et  al.  (1989). 
Similarly,  the  difference  of  38  cm  ±22  cm  in  Figure  4.7  and  20  cm  ±  22  cm  in  Figure  4.5 
between  the  German  height  datum  (NN)  and  French  height  system  (IGN69)  compares  well  with 
the  estimate  of  50  cm  given  by  Boucher  (1993,  personal  communication)  based  on  the 
comparison  of  heights  at  stations  ^ong  the  French-German  border.  Though  there  is  a  conflict  in 
the  Australian  Height  datum  (AHD71)  which  is  shown  as  72  cm  below  the  German  height 
system  (NN)  in  Rapp  (1993)  and  it  is  74  cm  above  the  German  height  datum  in  Fipre  4.5,  the 
difference  between  the  German  height  datum  (NN)  and  the  NA  VD88  datum  shown  in  Figure  4.5 
agrees  well  with  the  estimates  given  by  Rapp  (1993)  to  ±4  cm. 

Also  the  difference  of  32  cm  (or  33  cm,  depending  on  the  solution)  between  the  defined 
Global  Vertical  Datum  surface  and  the  adopted  reference  ellipsoid  suggests  that  the  exact  value 
of  the  semi-major  axis  'a'  should  be  about  6378136.62  m  or  6378136.63  m  which  agrees  well 
with  the  recent  estimates  of  'a'  from  TOPEX/Poseidon  data  (Rapp,  1993,  personal 
communication). 
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CHAPTERS 


GLOBAL  VERTICAL  DATUM  DEFINITION  -  A  SIMPLIFIED  APPROACH 


In  Chapter  3  we  have  considered  an  ideal  solution  for  defining  a  Global  Vertical  Datum 
(GVD),  assuming  that  the  required  data  are  available  everywhere  and  idso  to  required  precision  to 
achieve  an  accuracy  of  ±10  cm  in  its  definition.  In  fact  reliable  and  accurate  gravimetric  ^ta  are 
available  only  in  small  portions  of  the  world  like  in  the  United  States,  western  Eurc^,  etc.  and 
may  i^uire  substantial  effort  and  coordination  in  the  p^m  of  an  International  organization  to  come 
up  with  implementation  of  GVD  that  can  meet  the  scientific  and  practical  ne^  of  the  geodetic 
community.  A  simplified  approach  in  realizing  a  GVD  has  been  discussed  in  Rapp  and 
Balasubrainania  (1992, Sections  3  and  7).  As  point^  out  in  the  report  (ibid.)  the  realized  Vertical 
Datum  will  not  have  die  sophistication  of  the  ideal  GVD  but  could  meet  the  various  practical  needs 
when  the  highest  accuracy  is  not  required.  Theo^  behind  such  a  simplified  approach, 
development  of  Height  Bias  models,  and  certain  numerical  tests  are  presented  in  this  chapter  and 
the  discussions  closely  follow  the  concepts  explained  in  Rai^  and  Balasubramania  (ibid.). 

5.1  Basic  Concents 

This  api^ach  is  mainly  dependent  on  the  widely  available  space  geodetic  techniques  such 
as  Global  Positioning  System  (GPS)  and  DORIS  tracking  netwo^  which  are  becoming,  and 
promise  to  remain  for  sometime,  some  of  the  most  important  and  readily  accessible  geodetic 
measurement  systems,  in  particular  GPS.  The  concept  underl^g  this  approach  is  to  let  rnUP  92 
or  any  other  geocentric  reference  fiame  realized  by  the  International  Earth  Rotation  Service  (lERS) 
define  a  fundamental  reference  system.  Knowing  the  transformation  parameters  that  will  allow  the 
conversion  fiom  the  individual  coordinate  system  that  is  associated  with  the  tracking  type  used, 
such  as  GPS,  DORIS,  etc.  to  the  ITRF  system,  tire  rectangular  coordinates  of  the  staticms  will  be 
known  in  a  center-of-mass,  properly  scaled  cocndinate  system.  Now  assunung  an  ideal  reference 
ellipsoid,  after  defiiung  its  ^uattreial  radius  and  flattening  (refer  to  equation  4.1),  the  ellipsoidal 
height  h  at  the  stations  can  Ire  computed  using  standard  techniques.  The  ellipsoidal  heights  or 
height  differences  are  then  used  in  c  snjuncticm  with  the  getnd  u:>dula!ion  informs  ion  to  'compute 
orthometric  heights  using  equation  2.1.  The  Okti  ^metric  heights  thus  computed  refer  to  a 
ccirtceptua!  surface  the  gr oid)  which  is  not  physically  realizajle  ^id  also  net  associated  with  any 
specific  nrean  sea  level.  The  geoid  undulation  information  required  at  the  stations  can  either  be 
obtained  from  regional  geoid  height  models  like,  GEOID93  (in  USA),  GSD91  (in  Canada),  NKG- 
89  (in  Scandinavian  countries),  etc.  if  available,  or  can  be  computed  through  Stokes' 
technique/LSC  techrtique  by  combining  potential  coefficient  information  from  high  degree  (360) 
model  with  surface  gravity  data  (discuss^  in  Chapter  3).  In  areas  where  either  of  Jre  above  two 
tnetiiods  is  not  feasible,  the  global  geopotential  models  available  today  that  are  complete  to  degree 
360  can  be  used.  One  such  nxxlrl  giving  moderately  good  results  is  tire  OSU91 A  niodel  described 
by  Rapp,  Wang,  Pavlis  (1991).  The  use  and  accuracy  of  the  OSU91A  nxidel  for  geoid 
det^mination  are  described  by  Rapp  (1992)  and  Rapp  and  Balasubramania  (1992). 

The  limitations  of  this  tirethod  relate  to  the  precision  in  which  the  ellipsoidal  heights  or 
height  differences  can  be  determined  in  geocentric  reference  frames  with  GPS  arid  DORIS  tracking 
networks  and  ttiso  on  the  accuracy  of  available  geoidal  information.  Widely  available  GPS 
observations  (height  determination  specifically),  precise  undulation  values  and  a  Height  Bias  motfel 
are  essential  for  implementing  this  approach. 
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5.2  Development  of  a  Height  Bias  Model 


Since  the  regional  vertical  datums  that  are  realized  and  used  today  do  not  have  the  ideal 
geoid  as  reference  surface,  the  oithometric  heights  computed  on  a  regional  vertical  datum  will  be 
different  from  the  (nthometric  height  that  is  computed  with  respect  to  an  ideal  geoid.  The  difference 
may  simply  be  a  constant  or  it  may  be  a  difference  driven  by  a  number  of  factors  including 
distortions  in  a  local  vertical  network  (Rapp,  1992).  Let  Hj  be  the  orthometric  height  computed 
using  the  relation. 

Hi  =  h-N  (5.1) 

where  h  is  the  ellipsoidal  height  of  the  station  measured  along  the  ellipsoidal  ntnmal  passing 
through  the  station  and  N  is  the  geoidal  undulation  which  is  the  separation  between  the  defined 
ellipsoid  and  the  ideal  geoid  Let  be  the  orthometric  height  computed  in  a  local  vertical  datum. 
Assume  a  linear  relation  between  Hq  and  as  follows: 


Hd  =  Hi  +  c(<j>.X) 


(5.2) 


where  c(<(>,X)  can  be  considered  as  a  correction  term  or  height  bias  term,  dependent  on  position  in 
the  datum  to  correct  tnthometric  heights  referred  to  an  ideal  reference  surface  to  the  regional  vertical 
reference  datum.  If  sufficient  points  are  available,  the  value  of  'c'  can  be  mapped  and  represented 
in  some  functional  form  depending  cn  the  size  of  the  area  and  the  behavior  of  'c'.  Once  this  bias 
function  is  determined  for  a  region,  the  ideal  orthometric  height  computed  based  on  discussion  in 
the  previous  section  can  be  cemverted  into  the  orthometric  height  in  the  local  datum. 

Development  of  a  Height  Bias  model  or  setting  procedures  f(x  mapping  the  Bias  function 
'c'  in  some  functional  form  requires  GPS  observations  (height  determination  speciBcally),  precise 
undulation  data  and  consistent  orthometric  heights  in  a  local  vertical  datum.  Currently  such 
informatics  is  available  in  the  United  States,  C^ada  and  in  some  pam  of  Eurc^.  On  our  r^uest, 
the  Geodetic  Survey  of  Canada  (A.  MainviUe,  personal  communication,  April  1992)  provided  the 
GPS  data  for  four  traverses:  South  Alberta  (106  stations),  North  Alberta  (51  stations).  Central 
Alberta  (52  stations),  and  the  Great  Slave  Lake  Area  (91  stations).  In  addition,  the  orthometric 
heights  for  all  the  GPS  stations  in  the  (Canadian  Geodetic  Vertical  Datum-28  ((XiVD-28)  and  the 
high  resolution  geoid  tiKxlel  GSD91  were  also  obtained.  For  the  United  States,  the  National 
Gecxietic  Survey  (D.Milben,  personal  communication,  Oct  1992)  provided  the  data  for  GPS 
traverses  in  Florida  (52  stations),  Virginia  GIOS  (62  stations),  and  in  Oregon  (44  stations).  They 
also  supplied  the  or^ometric  heights  of  all  the  traverse  stations  in  National  Geodetic  Vertietd 
Datum  29  (NGVD-29)  and  the  high  resolution  Geoid  Height  Models  GEOID90  and  the  latest 
GEOID93.  GPS  traverse  data  for  49  stations  in  the  Tennessee  area  along  with  their  orthometric 
heights  received  from  Zeigler,  State  of  Tennessee,  Department  of  Transportation,  Nashville, 
Tennessee  were  also  used  for  developing  Height  Bias  model  in  the  Tennessee  area.  The 
VERTCON  (Version  1.00,1992)  software  (mentioned  in  Chapter  4)  received  from  Venical 
Network  Branch,  NGS  was  used  to  convert  NGVD-29  normal  orthometric  heights  to  NAVD88 
Helmert  orthometric  heights. 

Based  on  the  data  received,  the  following  procedure  was  followed  to  develop  a  suitable  Height 
Bias  model: 

(a)  Since  the  GPS  static  ii  eexmlinates  in  the  Canadian  GPS  traverses  as  well  as  Tennessee  GPS 
traverse  were  given  in  the  WGS84  ccwrdinate  system,  they  were  transformed  to  ITRF90 
geocentric  coordinate  system,  using  the  transformation  parameters  given  in  Boucher  and 
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Altamimi  (1991),  to  make  the  origin  of  the  coordinate  system  compatible  with  die  geopotoidal 
model  used  In  the  case  of  GPS  traverses  received  from  NGS,  the  station  coordinates  were 
already  in  the  ITRF90  geocentric  system  as  reported  by  Milben  (1991b)  and  no  transformation 
was  required.The  geoidal  undulations  were  tlwn  interpolated  using  the  cmresponding  geoid 
height  models. 

(b)  The  oithoroetric  height  bias  'c'  was  then  computed  using  the  relation 

c(<|>A)  =  (hGPS-NGEOID  model) -HlocAL  DATUM  (5-3) 

Exaii^les  of  mihometric  height  bias  tables  are  given  in  Table  S.l  for  die  Oregtxi  area  in  the 
USA  and  in  Table  5.2  for  the  Great  Slave  Lake  Area  in  Canada.  For  the  U.S  travo^s,  both  the 
currently  used  NGVD-29  heights  as  well  as  the  newly  established  NAVD88  heights  were  used  in 
developing  the  Height  Bias  functions.  Statistical  aiialysis  of  the  height  bias  values  for  U.S  and 
Canadian  traverses  are  given  in  Table  5.3. 

Table  5.3  Statistics  of  the  Difference  in  Orthometiic  Heights  Derived  From  GPS/GSD91  and 
CGVD28  Heights  in  Case  of  Canadian  Traverses  and  From  GPS/GEOID93  and 
NAVD88  Heists  for  GPS  Traverses  in  the  United  States.  Units  are  in  cm. 
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Table  S.  1  Computation  of  OrthcHnetric  Height  Bias  Function  (c) 
Oregon.  United  States.  Heights  are  in  meters. 

,  in  Oregon  GPS  Traverse  Area, 

Lat 

(deg) 

Ell.ht 

(h) 

Ortho. ht 
(H) 

(HGVD29} 

c  -  {h-N)-H 

Ortho. ht  Geoid.ht  Ht.bias  value 
(H)  (N)  (tlAVD88)  (N6VD29) 

(NAVD88)  {GB0I093)  (units  in  ca) 

45.47 

239.26 

557.101 

577.160 

578.180 

-20.573 

-50.6 

51.4 

42.22 

237.29 

533.994 

556.686 

557.736 

-22.988 

-75.4 

29.6 

44.83 

235.93 

-6.483 

16.121 

17.121 

-22.817 

-78.7 

21.3 

43.98 

242.76 

665.492 

682.364 

683.344 

-17.630 

-22.2 

75.8 

43.13 

238.20 

1392.601 

1412.875 

1414.165 

-20.648 

-91.6 

37.4 

43.81 

239.40 

1394.157 

1412.739 

1413.949 

-19.044 

-74.8 

46.2 

45.52 

237.01 

36.439 

58.039 

59.089 

-21.098 

-74.2 

30.8 

46.11 

236.79 

-16.333 

3.562 

4.542 

-20.297 

-57.8 

40.2 

45.30 

242.18 

845.869 

862.066 

863.176 

-17.049 

-25.8 

85.2 

44.39 

236.71 

59.025 

80.570 

81.590 

-21.900 

-66.5 

35.5 

45.24 

239.82 

863.124 

882.423 

883.513 

-19.975 

-41.4 

67.6 

45.00 

239.78 

784.779 

803.820 

804.930 

-19.616 

-53.5 

57.5 

43.74 

241.92 

882.872 

899.963 

901.063 

-17.773 

-41.3 

68.7 

42.09 

241.43 

1255.053 

1273.403 

1274.553 

-19.139 

-36.1 

78.9 

45.86 

237.20 

1.696 

22.843 

23.783 

-21.548 

-53.9 

40.1 

42.98 

242.95 

1323.322 

1339.038 

1340.108 

-16.443 

-33.8 

73.2 

45.35 

242.75 

1246.949 

1261.640 

1262.940 

-15.771 

-22.0 

108.0 

42.32 

237.12 

400.721 

423.829 

424.849 

-23.335 

-79.3 

22.7 

45.61 

238.81 

17.243 

37.637 

36.637 

-20.810 

-58.4 

41.6 

43.23 

242.74 

1971.917 

1988.030 

1989.070 

-16.749 

-40.4 

63.6 

45.13 

236.03 

12.240 

33.989 

35.029 

-22.019 

-77.0 

27.0 

44.93 

236.77 

41.651 

62.907 

63.927 

-21.539 

-73.7 

28.3 

46.30 

235.92 

53.461 

76.690 

77.720 

-23.670 

-58.9 

44.1 

43.12 

235.82 

3.745 

27.942 

29.012 

-24.483 

-78.4 

28.6 

42.70 

239.45 

1334.608 

1354.065 

1355.265 

-19.834 

-82.3 

37.7 

44.17 

235.88 

-15.986 

7.050 

8.100 

-23.219 

-86.7 

18.3 

45.71 

238.45 

81.190 

101.612 

102.632 

-20.773 

-66.9 

35.1 

44.44 

241.30 

1099.198 

1115.895 

1117.105 

-17.553 

-35.4 

85.6 

43.70 

235.90 

-21.098 

2.492 

3.562 

-23.925 

-73.5 

33.5 

45.80 

236.54 

416.104 

435.452 

436.522 

-19.856 

-56.2 

50.8 

45.47 

236.16 

-15.598 

5.530 

6.580 

-21.365 

-81.3 

23.7 

43.25 

236.64 

140.826 

163.290 

164.350 

-22.788 

-73.6 

32.4 

44.43 

238.05 

1116.756 

1136.953 

1138.153 

-20.494 

-90.3 

29.7 

45.48 

236.30 

17.992 

38.079 

39.139 

-20.419 

-72.8 

33.2 

44.31 

238.44 

947.275 

967.450 

968.600 

-20.483 

-84.2 

30.8 

43.59 

240.93 

1256.144 

1273.921 

1275.111 

-18.452 

-51.5 

67.5 

44.55 

239.96 

1301.683 

1319.711 

1320.941 

-18.798 

-46.0 

77.0 

44.15 

237.43 

221.158 

242.724 

243.854 

-21.958 

-73.8 

39.2 

42.19 

239.64 

1425.684 

1445.488 

1446.678 

-20.060 

-93.4 

25.6 

43.01 

239.23 

1328.545 

1347.868 

1349.088 

-19.742 

-80.1 

41.9 

43.01 

236.30 

302.623 

325.620 

326.680 

-23.320 

-73.7 

32.3 

43.64 

236.42 

21.257 

43.451 

44.541 

-22.538 

-74.6 

34.4 

42.47 

235.84 

1126.284 

1151.520 

1152.590 

-25.296 

-101.0 

6.0 

42.01 

242.28 

1343.711 

1361.257 

1362.397 

-18.209 

-47.7 

66.3 
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Table  S.2  Computation  of  Orthometric  Height  Bias  Function  in  Great  Slave  Lake  Area.  Canada. 
Heights  are  in  meters. 


L«t 

(deg) 

Lon  snip. he  goold.ht  (h-N)  Ortho. ht  Ht.orror 

(dtg)  (h)  (M)  (H)  (in  ea) 

(froB  0SD91)  (in  CGVD28) 

62.475 

245.559 

183.350 

-25.905 

209.255 

209.229 

2.5 

60.038 

243.117 

279.969 

-19.840 

299.809 

299.902 

-9.2 

60.426 

243.647 

263.951 

-21.571 

285.522 

285.730 

-20.7 

60.800 

243.408 

244.827 

-21.935  , 

266.762 

267.028 

-26.5 

60.709 

244.995 

216.481 

-24.809 

241.290 

241.269 

2.0 

62.479 

245.273 

155.184 

-25.623 

180.807 

180.706 

10.1 

62.544 

245.010 

141.715 

-25.083 

166.798 

166.700 

9.8 

62.593 

244.813 

142.511 

-24.582 

167.093 

167.061 

3.2 

62.660 

244.705 

152.975 

-24.274 

177.249 

177.221 

2.8 

62.686 

244.535 

148.808 

-23.980 

172.788 

172.753 

3.4 

62.709 

244.411 

147.611 

-23.776 

171.387 

171.339 

4.8 

62.784 

244 . 007 

140.531 

-22.935 

163.466 

163.455 

1.0 

62.641 

243.743 

210.678 

-22.284 

232.962 

232.996 

-3.3 

62.556 

243.591 

233.449 

-22.033 

255.482 

255.507 

-2.5 

62.382 

243.503 

250.968 

-21.915 

272.883 

272.957 

-7.3 

62.302 

243.571 

221.520 

-22.063 

243.583 

243.699 

-11.6 

62.204 

243.683 

206.345 

-22.334 

228.679 

228.722 

-4.3 

62.109 

243.704 

191.012 

-22.422 

213.434 

213.524 

-9.0 

62.024 

243.687 

188.108 

-22.451 

210. 5S9 

210.650 

-9.1 

61.979 

243.563 

196.285 

-22.232 

218.517 

218.608 

-9.1 

61.837 

243.329 

203.860 

-21.826 

225.686 

225.825 

-13.8 

61.781 

243.228 

201.425 

-21.666 

223.091 

223.148 

-5.7 

61.711 

243.103 

197.065 

-21.465 

218.530 

218.533 

-0.3 

61.677 

242.975 

190.166 

-21.261 

211.427 

211.571 

-14.4 

61.601 

242.859 

175.645 

-21.128 

196.773 

196.921 

-14.7 

61.514 

242.754 

154.450 

-21.025 

175.475 

175.646 

-17.1 

61.428 

242.601 

138.736 

-20.823 

159.559 

159.752 

-19.2 

61.366 

242.499 

136.108 

-20.732 

156.840 

157.021 

-18.0 

61.095 

242.500 

174.119 

-20.585 

194.704 

194.958 

-25.4 

60.870 

243.269 

206.432 

-21.749 

228.181 

228.496 

-31.4 

60.717 

243.528 

231.687 

-22.030 

253.717 

253.948 

-23.0 

60.629 

243.659 

240.712 

-22.095 

262.807 

262.961 

-15.3 

60.784 

245.342 

211.651 

-25.497 

237.148 

237.128 

2.0 

60.611 

245.510 

245.756 

-25.469 

271.225 

271.277 

-5.2 

60.536 

245.610 

222.953 

-25.550 

248.503 

248.527 

-2.4 

60.407 

245.733 

242.665 

-25.609 

268.274 

268.213 

6.1 

60.228 

246.153 

247.940 

-25.958 

273.898 

273.778 

11.9 

60.174 

246.295 

239.781 

-26.070 

265.851 

265.705 

14.6 

60.149 

246.444 

241.794 

-26.275 

268.069 

267.889 

18.0 

60.120 

246.619 

232.888 

-26.527 

259.415 

259.224 

19.0 

60.035 

246.872 

226.431 

-26.810 

253.241 

253.042 

19.9 

60.026 

247.027 

219.796 

-27.054 

246.850 

246.597 

25.3 

60.014 

247.753 

161.315 

-28.010 

189.325 

188.848 

47.7 

60.013 

247.955 

172.707 

-28.241 

200.948 

200.560 

38.7 

60.005 

248.046 

178.882 

-28.335 

207.217 

206.868 

34.8 

59.999 

248.162 

183.599 

-28.439 

212.038 

211.711 

32.7 

62.458 

245.413 

165.830 

-25.813 

191.643 

191.579 

6.3 
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Tabic  5.2  (cont.) 


62.519 

245.101 

160.654 

-25.285 

185.939 

185.810 

12.9 

62.759 

244.193 

158.507 

-23.318 

181.825 

181.830 

-0.4 

60.456 

245.675 

229.981 

-25.583 

255.564 

255.514 

4.9 

60.283 

245.955 

245.985 

-25.714 

271.699 

271.624 

7.4 

61.185 

246.301 

132.292 

-27.566 

159.858 

159.663 

19.4 

61.313 

242.400 

138.238 

-20.580 

158.818 

158.955 

-13.6 

60.106 

243.242 

269.264 

-20.182 

289.446 

289.596 

-14.9 

60.173 

243.308 

267.289 

-20.420 

287.709 

287.873 

-16.3 

60.263 

243.424 

269.261 

-20.821 

290.082 

290.243 

-16.0 

60.343 

243.555 

269.673 

-21.223 

290.896 

291.060 

-16.3 

60.506 

243.747 

251.176 

-21.990 

273.166 

273.324 

-15.8 

60.625 

243.939 

190.459 

-22.742 

213.201 

213.289 

-8.7 

60.684 

244.067 

170.230 

-23.100 

193.330 

193.446 

-11.6 

60.739 

244.482 

165.580 

-23.981 

189.561 

189.673 

-11.1 

60.725 

244.645 

167.320 

-24.253 

191.573 

191.611 

-3.8 

60.740 

245.179 

204.241 

-25.163 

229.404 

229.384 

2.0 

60.804 

245.465 

206.889 

-25.756 

232.645 

232.681 

-3.5 

61.039 

246.312 

135.807 

-27.525 

163.332 

163.270 

6.1 

60.838 

244.220 

138.667 

-23.710 

162.377 

162.663 

-28.6 

60.000 

247.395 

197.123 

-27.534 

224.657 

224.418 

23.8 

60.047 

247.190 

203.371 

-27.345 

230.716 

230.463 

25.2 

60.046 

246.745 

242.739 

-26.620 

269.359 

269.141 

21.7 

60.329 

245.851 

241.294 

-25.623 

266.917 

266.865 

5.2 

60.670 

245.404 

232.053 

-25.382 

257.435 

257.380 

5.5 

60.698 

245.256 

230.009 

-25.197 

255.206 

255.157 

4.8 

60.825 

245.629 

190.050 

-26.062 

216.112 

216.110 

0.2 

60.865 

245,759 

179.217 

-26.381 

205.598 

205.594 

0.4 

60.928 

245.883 

165.257 

-26.669 

191.926 

191.923 

0.2 

60.986 

246.180 

139.620 

-27.235 

166.855 

166.851 

0.4 

61.132 

246.370 

131.929 

-27.653 

159.582 

159.459 

12.3 

60.965 

246.045 

137.190 

-27.013 

164.203 

164.180 

2.346 

60.731 

244.844 

178.432 

-24.599 

203.031 

203.107 

-7.581 

60.742 

244.260 

151.468 

-23.571 

175.039 

175.221 

-18.233 

60.766 

244.148 

148.883 

-23.434 

172.317 

172.526 

-20.886 

60.559 

243.860 

236.930 

-22.382 

259.312 

259.488 

-17.575 

60.001 

243.018 

276.084 

-19.586 

295.670 

295.814 

-14.376 

60.985 

242.753 

171.685 

-20.924 

192.609 

192.884 

-27.462 

60.931 

243.082 

226.397 

-21.475 

247.872 

248.167 

-29.489 

60.936 

242.933 

214.140 

-21.195 

235.335 

235.559 

-22.448 

61.048 

242.599 

187.730 

-20.697 

208.427 

208.678 

-25.071 

61.173 

242.469 

152.332 

-20.634 

172.966 

173.093 

-12.733 

61.255 

242.473 

134.933 

-20.709 

155.642 

155.804 

-16.164 

61.731 

243.139 

208.827 

-21.516 

230.343 

230.492 

-14.921 

60.654 

243.634 

243.395 

-22.101 

265.-496 

265.708 

-21.137 
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Though  several  procedures  exist  to  map  the  bias  ^inction  'c',  the  following  three  are  the  generally 
used  techniques  for  interpolating  a  function  value  at  any  point,  given  the  function  values  at  certain 
known  points  in  that  area.  A  few  other  techniques  are  also  discussed  in  Shresdia  et  al  (1993). 

(a)  Global  surface  fitting  using  deterministic  polynomial  functions  :  This  technique  provides  a 
simple  way  of  modelling  the  bias  function  value.  Different  types  of  polynomial  functions  can 
be  used  to  model  the  Height  Bias  function  'c',  for  example  -  bilinear  polynomials,  quadratic 
polynomials,  cubic  polynomials,  etc.  Though  it  is  a  very  simple  technique,  these  polynomial 
functions  provide  unpr^ctable  errors  in  the  interpolated  values  away  from  the  known  station 
points,  m^ng  it  unsuitable  for  our  purpose. 

(b)  Minimum  curvature  fitting  using  spline  functions  :  Cubic  spline  functions  provide  a 
continuous  and  smooth  composite  function  fitting  the  known  data  more  accurately  without  the 
undue  oscillations  caused  by  polynomial  functions  away  from  the  known  stations.  Due  to  the 
versatility  of  this  technique,  it  has  been  implennented  in  the  NADCON  software  of  the  NGS 
(Dewherst,  1990)  for  the  transformation  of  positional  data  firom  the  NAD27  to  NAD83 
coordinate  system.  The  main  disadvantage  of  this  technique  is  that  the  accuracy  of  the 
predicted  value  of  the  function  cannot  be  estimated  in  a  simple  way.  Though  the  height  bias 
plots  made  from  this  method  for  (Canadian  GPS  traverse  areas  showed  a  goc^  agreement  with 
plots  made  using  least  ^uares  collocation  technique,  it  was  not  chosen  to  mo^l  the  height 
bias  function  for  the  main  reason  that  the  accuracy  of  the  predicted  bias  function  value  couldn't 
be  determined  in  a  simpler  way  as  in  the  case  of  the  LSC  technique. 

(c)  Height  Bias  modelling  using  least  squares  collocation  technique  :  As  explained  in  Chapter  3, 
this  technique,  in  addition  to  providing  the  Best  Linear  Prediction  of  height  bias  values,  also 
gives  the  accuracy  of  the  pre^cted  vdue  as  given  in  equation  (3.76)  though  sometimes  too 
optimistic.  This  technique  was  chosen  to  develop  Height  Bias  model  functions  in  the  test 
areas  where  GPS  traverse  data  and  geoid  undulation  data  were  made  available. 

5.3  Numerical  Results 

Assuming  a  correlation  length  of  40  km  and  data  noise  of  ±5  cm,  with  a  second-order 
Markov  model  covariance  function  fitted  to  the  kix>wn  height  bias  data  (after  removing  the  mean)  at 
the  GPS  traverse  stations,  the  Height  Bias  value  along  with  its  accuracy  was  predicted  for  any 
other  location  in  the  area  using  least  squares  collocation  technique.  (Contour  plots  showing  the 
predicted  Height  Bias  function  in  the  area  and  its  accuracy,  one  each  for  U.S  traverse  (in  Oregon) 
and  Canadian  traverse  (in  South  Alberta)  are  shown  in  Figures  S.l,  5.2  and  5.4, 5.5,  respectiv 
with  a  contour  interval  of  5  cm.  The  perspective  view  of  the  Height  Bias  functions  are  also  shown 
as  5.3  and  5.6  for  those  areas.  To  compare  the  Height  Bias  function  values  with  different  geoid 
height  models  as  well  as  height  systems  (NAVD88  and  NGVD-29)  two  other  plots  are  also 
included  as  Figures  5.7  and  5.8.  An  orthometric  height  bias  function  plot  in  case  of  the 
Scandinavian  GPS  traverse  computed  based  on  normal  heights  of  the  stations  derived  from  UELN- 
73  geopotential  numbers  is  also  given  as  Figure  5.9. 

g.4  Discussion  of  Results 

From  Figures  5.1  and  5.4,  we  see  the  sulvantage  of  developing  such  Height  Bias  models. 
By  knowing  the  orthometiic  height  with  respect  to  an  ideal  geoid  using  Geoid/GPS  observations 
using  equation  (5.1),  its  orthcmetric  height  in  a  local  vertical  datum  can  be  computed  using  the 
interpolated  bias  function  'c'.  Also  the  plot  showing  the  accuracy  of  predicted  'c'  values  serves  as 
a  reliability  diagram  intimating  the  user  with  the  accuracy  of  the  predicted  'c'  value  and  hence  the 
orthometric  height  in  the  local  vertical  datum.  From  Figures  5.2  and  5.5,  we  see  the  accuracy  of 
predicted  height  bias  function  is  of  the  order  of  ±20  cm  if  we  have  well  distributed  GPS  stations  in 
the  area.  Considering  Figure  5.9,  one  can  see  the  Height  Bias  function  modelled  for  the 
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Scandinavian  GPS  traverse,  which  runs  through  four  different  countries,  namely,  Germany, 
Denmark,  Sweden,  and  Norway  cannot  find  its  intended  use  for  the  following  reasons; 

•  Height  Bias  function  which  was  computed  based  on  normal  heights  of  the  stations  derived 
from  UELN-73  geopotential  numbers  is  not  compatible  with  the  National  height  systems  used 
in  different  countries  involved  in  the  GPS  traverse. 

•  The  GPS  traverse  is  almost  a  straight  line  traverse  running  from  South-West  to  North-East 
direction  and  covering  no  area  as  such. 

One  can  conclude  that  this  plot  indicates  how  a  GPS  traverse  should  not  be  used  if  we  want  to 
develop  a  Height  Bias  model  for  an  area  with  available  undulation  data.  Comparing  the  contour 
plots  in  Figures  S.l  and  5.7,  the  clear  signature  in  the  eastern  part  of  Oregon  due  to  the  large  bias 
terms  in  using  NGVD-29  heights,  can  be  seen  to  be  minimized  when  we  use  NA  VD88  heights. 
But  we  see  the  large  bias  terms  in  the  western  part  of  Oregon  exist  for  both  the  height  systems.  The 
possible  reason  for  this  may  be  the  error  in  leveling  data  caused  by  the  network  design  in 
mountainous  areas  or  may  be  due  to  the  sea  level  slope  observed  along  the  Pacific  Coast.  As 
pointed  out  by  Zilkoski  et  al  (1992),  there  are  long  north-south  leveling  lines  through  the  valleys, 
and  only  a  few  east-west  lines  that  cross  the  mountains  to  create  loops  to  help  control  and  check 
accumulation  of  errors. 

5.5  Summary  and  Conclusions 

1 .  The  deHnition  of  a  Global  Vertical  Datum  (GVD)  using  the  procedures  described  in  this 
Chapter  eliminate  a  direct  reference  to  mean  sea  level.  This  means  that  countries  that  do  not 
have  access  to  the  oceans  can  also  use  such  a  GVD  if  precise  geocentric  coordinates  for  one  or 
more  points  can  be  determined  in  the  country  (Rapp,  1983). 

2.  The  accuracy  limitation  of  this  approach  relates  to  the  precision  in  which  the  ellipsoidal  heights 
or  height  differences  can  determined  in  the  geocentric  reference  frame  with  GPS  or  any  other 
space  geodetic  techniques  and  also  on  the  accuracy  of  geoidal  information  available. 

3.  For  getting  the  undulation  information  we  have  assumed  that  the  terrestrial  gravity  has  had  a 
common  teight  system  in  the  calculation  of  gravity  anontalies.  This  assun^tion  is  cleariy  false 
and  hinders  the  fundamental  goal  at  the  ±10  cm  level  of  accuracy  (Rapp  and 
Balasubramania,  1 992). 

4.  The  simplified  approach  discussed  in  this  Chapter  does  provide  a  means  for  determining  the 
Global  Vertical  Datum  since  the  geoid  forms  the  reference  surface.  Unfortunately,  we  lack 
sufficiently  precise  knowledge  of  the  geoid  undulations  all  over  the  globe  to  make  the 
procedure  universally  applicable. 

5.  The  least  squares  collocation  (LSC)  technique  described  in  Chapter  3,  provides  an  accurate 
procedure  for  modelling  the  Height  Bias  function  and  also  to  compute  the  accuracy  of  such 
Height  Bias  function. 

6.  The  modelled  Height  Bias  function  refers  only  to  the  particular  geoid  model  and  local  height 
system  used  in  its  ^velopment. 

7 .  Well  distributed  control  points  in  the  local  height  datum  help  in  developing  a  meaningful  Height 
Bias  function  that  would  enable  the  orthometric  heights  computed  with  respect  to  the  geoid  be 
converted  to  the  specific  vertical  datum  system. 
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HEIGHT  BIAS  FUNCTION  PLOT 


Contour  Interval  5  cm 
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I.ONCITUOE 


Figure  5.1  Height  Bias  Function  'c*  Interpolated  Using  Least  Squares  CoUocadon 
Technique  with  Oregon  GPS  Traverse  Data  Set  (GEOID93/NAVD88) 
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LATITUDE 


ACCURACY  ESTIMATE  OF  BIAS  FUNCTION  PLOT 
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Figure  5.2  Accuracy  Plot  of  Height  Bias  Function  'c'  Inteipolated  Using  Least  Squares 
Collocation  Technique  with  Oregon  GPS  Traverse  Data  Set 
(GEOID93/NAVD88) 
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JiEIGHT  BIAS  FUNCTION  -  A  PERSPECTIVE  VIEW 


Figure  S.3  Perspective  View  of  the  Height  Bias  Function  'c'  Interpolated  Using  Least 
Squares  Collocation  Technique  with  Oregon  GPS  Traverse  Data  Set 
(GEOID93/NAVD88) 
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HEIGHT  BIAS  FUNCTION  PLOT 
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Figure  5.4  Height  Bias  Function  ’c’  Interpolated  Using  Least  Squares  Collocation 
Technique  with  South  Alberta  (Canada)  GPS  Traverse  Data  Set 
(GSD91/CGVD28) 


88 


LATITUDE 


LfiTlTJOE 


ACCURACY  OF  HEIGHT  BIAS  FUNCTION  PLOT 
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Figure  S.S  Accuracy  Plot  of  Height  Bias  Function  'c*  Interpolated  Using  Least  Squares 
Collocation  Technique  with  South  Alberta  (Canada)  GPS  Traverse  Data  Set 
(GSD91/CGVD28) 
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LRTITJDE 


HEIGHT  BIAS  FUNCTION  -  A  PERSPECTIVE  VIEW 


CIKECTION  OF  VIEW 


Fieure  5.6  Perspective  View  of  the  Height  Bias  Function  'c'  Interpolated  Using  Least 
Squares  Collocation  Technique  with  South  Alberta  (Canada)  GPS  Traverse 
Data  Set  (GSD91/CGVD28) 
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HEIGHT  BIAS  FUNCTION  PLOT 
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Figure  5.7  Height  Bias  Function  'c'  Interpolated  Using  Least  Squares  Collocation 
Technique  widi  OregtHi  GTC  Traverse  Data  Set  (GEOID93/NGVD29) 
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Figure  5.8 


Height  Bias  Function  'c'  Interpolated  Using  Least  Squares  Collocation 
Tec^que  with  Oregcm  GPS  Traverse  Data  Set  (GEOID90/NGVD29) 
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Hgure  5.9  Height  Bias  Function  Plot  using  Scandinavian  GPS  Traverse  Data/UELN 
Heights 
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CHAFreR6 


SUMMARY  AND  CONCLUSIONS 


Extending  the  concepts  prc^tosed  by  Rapp  and  Balasubramania  (1992)  for  fonnulating  the 
World  Height  System,  this  study  has  considered  two  different  approaches  for  defining  and 
realizing  a  Global  Vertical  Datum  (G  VD).  The  first  approach  is  an  idealistic  approach  in  which 
accurate  infonnation  from  various  regitxial  vertical  datums  are  brought  together  in  an  adjustment 
leading  to  the  definition  of  GVD  and  the  other  one  is  a  simplifi^  i4>proach  based  on  space 
positioning  and  geoid  undulation  information  with  geoid  as  reference  surface.  While  the  geoid  in 
the  latter  approach  has  no  direct  relationship  with  a  mean  sea  level,  the  regional  vertical  datum 
information  used  in  the  idealistic  approach  will  be  dependent  on  the  mean  sea  level  observatiois  at 
one  ot  more  tkle  gauges  for  their  de^tion  and  estabtishment. 

After  a  brief  review  of  die  different  heights  and  height  systems  practically  used  in  different 
parts  of  the  world  along  with  their  intercomparisons,  the  current  procedures  used  to  define  and 
establish  a  regional  vertical  datum  using  the  long  term  mean  sea  level  observations  at  one  or  more 
primary  tide  gauge  stations  were  presented.  It  was  em|diasized  here  that  since  die  mean  sea  levels 
computed  at  the  different  tide  gauges  do  not  coincide  with  the  geoid  due  to  the  existence  of  Sea 
Surface  Topography  (SST)  caus^  by  ocean  currents  and  various  other  factors,  the  resulting 
adjustment  of  leveling  networks  to  fit  diese'mean  sea  tevels  will  result  in  distortions  in  the  realized 
regional  vertical  datums.  Several  researchers  have  come  up  with  ways  and  means  to  remove  such 
distortions  including  attempts  to  model  the  SST  at  the  tide  gauges  as  an  absolute  value  and  also  as 
difference  between  adjacent  tide  gauges.  After  listing  the  efforts  so  far  made  by  numerous 
researchers,  a  mathematical  model  for  the  combined  adjustment  solution  to  define  an  id^  Regional 
Vertical  Datum  was  proposed.  This  prc^iosed  solution  was  based  on  the  assumption  that  realistic 
models  to  estimate  tlw  SST  values  at  the  tide  gaug^  are  available.  But  since  we  realize  that  perfect 
noodelling  of  SST  at  the  tide  gauges  is  not  possible,  at  least  in  the  near  future,  and  there  will  be 
inconsistencies  in  the  estimt^d  SST  values  at  the  tide  gauges,  adopting  the  robust  sequential 
procedure  put  fOTward  by  Schafhin  (1986, 1990)  was  recommended  for  re^zing  an  ideal  regional 
vertical  datum. 

Mathematical  modelling  required  to  relate  the  pritiuuy  observables  to  the  parameters  to  be 
established  was  studied  for  Afferent  data  scenarios.  Separate  modelling  procedures  were 
developed  for  setting  up  the  observatitm  equations  at  the  space  geodetic  stations,  depending  on  die 
type  of  fiee-air  anomaly  data  and  heights  used  in  a  regional  vertical  datum.  Special  emphasis  was 
also  placed  on  identifying  different  etm'  sources  affecting  the  terrestrial  gravity  anomalies  used  for 
defining  and  realizing  a  GVD  to  an  accuracy  of  ±10  cm.  The  effect  of  horizontal  datum 

inconsistencies  in  the  used  gravity  data  (Sg)])  and  the  corresponding  correction  needed  to  the 

computed  gravimetric  height  anomaly/  undulation  were  also  determined.  The  maximum  Sg^ 
correction  amounted  to  -0.097  mgal  at  the  Yarragadee  station  in  Australia  with  a  resulting 
undulation  coirectitx)  of  about  -2.3  cnt  Though  case  studies  have  shown  that  for  a  highly  accurate 
computation  of  height  anomaly  using  equation  (3.24),  consideratitm  of  Molodensky's  correction 
terms  up  to  g2  was  essential,  as  dense  grid  of  gravity  data  as  required  for  such  computations  may 
not  be  available  all  over  the  globe  in  the  near  future,  correction  terms  only  up  to  gi  were 
ctmsidered  sufficient  in  the  modelling  procedure. 
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A  simplified  approach  fm  defining  the  GVD  was  also  proposed  in  which  the  geoid  forms 
the  ideal  refovnce  surface.  Widely  available  GPS  observations  (height  determination  specifically), 
precise  undulation  values,  and  a  Height  Bias  model  were  consider^  essential  for  this  approach. 
Procedures  for  developing  a  Height  Bias  function  which  could  enable  the  orthometric  heights 
computed  with  respect  to  the  geoid  be  converted  to  the  specific  vertical  datum  system  were  also 
discussed.  Height  Bias  functions  developed  for  certain  areas  along  with  their  pr^cted  accuracy 
were  also  presented  in  the  form  of  contour  plots. 

Two  different  techniques,  namely  modified  Stokes'  technique  and  least  squares  collocation 
technique  (LSC),  were  consisted  for  computing  the  ^vimetric  height  anomaly/  undulation  at  the 
space  geodetic  stations  using  the  surface  ^vity  data  in  a  small  cap  around  the  computation  point 
and  potential  coefficient  informadon  from  a  global  geopotential  model.  To  avoid  dme  consuming 
and  unstable  numerical  inversion,  the  number  of  points  considered  for  predicting  gravimetric 
height  anomaly/undulation  using  LSC  technique  was  restricted  to  1600  points. 

For  numerical  analysis  OSU91 A  nrodel  (to  degree  360)  was  considered  as  reference  model. 
The  recently  made  available  JGM-2  (NASA  model)  from  degree  2  to  70,  and  augmented  by 
OSU91 A  niodel  coefficients  from  degree  71  to  360  was  also  used  to  estimate  the  remote  zone 
contribution  (beyond  the  cap  radius  in  which  surface  gravity  data  were  used)  for  the  computation 
of  gravimetric  height  anomaly/undulation  using  modified  Stokes'  technique.  The  height 
anomaly/undulations  computed  using  the  OSU91A  reference  model  differed  from  the  one 
computed  using  JGM-2  (augmented)  model  by  2  to  17  cm  at  the  space  geodetic  stations  in  the 
United  States,  4  to  25  cm  at  the  European  and  Scandinavian  stations  and  about  5  to  9  cm  at  the 
Australian  stations.  These  differences  are  consistent  with  the  accuracy  estimate  of  height 
anomalyAindulatitms  shown  in  Table  4.19  for  the  Modified  Stokes'  technique. 

A  first  iteration  attempt  to  estimate  the  parameters  defining  the  GVD  was  also  carried  out 
using  the  currently  available  data  in  six  regional  vertical  datums.  Data  used  in  the  numerical  studies 
had  ^e  following  shortcomings: 

•  Limited  gravity  data  were  available  in  the  Scandinavian  datum,  sufficient  enough  to  use  the 
gravity  data  in  a  cap  radius  of  0*.S  to  T.O  only,  though  capsize  of  2*  was  ctxisidered  optimal. 

•  Required  data  were  available  only  at  very  few  stations  in  the  various  regional  datums 
considered.  While  the  number  of  stations  us^  in  the  Australian  datum  was  two,  it  was  limited 
to  just  one  in  Germany,  France  and  England  regional  vertical  datums.  Also  the  gravity  data 
that  were  used  at  these  stations  was  not  corrected  for  terrain. 

•  In  the  United  States,  the  NAVD88  heights  used  did  not  result  from  direct  precise  leveling 
measurements  from  an  adjusted  Bench  Mark.  They  were  interpolated  using  the  VERTCON 
(Version  1.00)  software,  which  gives  the  corrections  to  convert  NGVD-29  normal  orthometric 
heights  to  NAVD88  Helmert  orthometric  heights,  only  up  to  cm  accuracy. 

The  following  conclusions  were  arrived  at  based  on  the  numerical  investigations  carried  out 
to  realize  a  GVD,  using  currently  available  data  with  the  mathematical  modelling  procedures 
developed  for  an  idealistic  approach; 

1 .  The  different  regional  vertical  surfaces  used  in  the  deflnition  of  GVD  show  the  same  pattern  of 
their  relative  positions  to  one  another,  in  spite  of  using  two  different  techniques,  such  as  the 
Modified  Stokes'  and  the  LSC  technique,  in  their  estimation. 

2 .  The  GVD  can  be  realized  to  an  accuracy  of  ±5  cm  with  the  accuracies  of  regional  vertical  datum 
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ccmnections  to  GVD  ranging  from  ±5  to  ±23  cm.  It  should  be  emphasized  that  the  error 
estimates  reported  here  ctmespond  to  a  'worst-case'  scenario,  where  ±7  mgal  was  the  accuracy 
assumed  for  the  terrestrial  gravity  data  available  at  the  Eurc^)ean  and  Australian  space  geodetic 
stations  and  the  variance-covariance  matrix  set  up  for  modified  Stokes'  technique  was  used  in 
the  computations. 

3.  The  results  (refer  to  Figure  4.5)  obtained  from  this  study  show  a  major  conflict  with  the  results 
reported  by  Rapp  (1993).  The  Australian  height  datum  AHD71  which  is  shown  as  74  cm 
above  the  German  height  datum  NN  in  Figure  4.5  of  this  report,  is  shown  as  72  cm  below  the 
German  height  datum  in  Rapp  (1993,  fig  1). 

4.  Other  than  the  discrepancy  stated  in  (3)  above,  the  estimated  separation  between  the  regional 
vertical  surfaces  compared  mostly  well  with  the  results  reported  by  various  geodesists  and 
oceanographers  based  on  their  detailed  regional  studies.  Some  of  the  compared  results  are  : 

•  Rossiter  ( 1967,  p  295)  based  on  his  analysis  of  annual  sea  level  variations  in  European  waters, 
has  repotted  that  the  average  value  of  the  mean  sea  level  in  the  northern  Mediterranean  (where 
the  Marseille  tide  gauge  is  located)  is  20  cm  below  the  N.A.P  (on  which  German  height  datum 
is  based).  TTiis  suggests  that  IGN  69  (France)  is  about  20  cm  below  NN  (Germany).  Boucher 
(personal  communication,  1993)  has  indicated  that  IGN  69  is  about  50  cm  below  NN 
(Germany).  The  result  from  this  study  shows  that  IGN  69  is  20  cm  ±22  cm  below  NN  datum. 

•  Also  based  on  oceanographic  studies,  Lisitzin  (1974)  shows  (refer  figure  30,  p  149,  ibid)  that 
the  mean  sea  level  at  southern  England  (where  the  Newlyn  tide  gauge  is  located)  is  atx>ut  20  cm 
below  the  mean  sea  level  in  southern  Fr^ce  suggesting  that  the  ODN  (England)  may  be  20  cm 
lower  than  the  IGN  69  (France).  Willis  et  al  (1989)  ba^  on  geodetic  observations  using  GPS 
techniques  have  reported  that  ODN  is  30cm  ±8cm  below  IGN  69.  The  result  from  this  study 
tdso  shows  that  ODN  is  33  cm  below  IGN  69. 

•  The  difference  between  German  height  datum  NN  and  NAVD88,  reported  as  76  cm  by  Rapp 
(1993),  agrees  well  with  72  cm  fiom  this  study. 

5.  The  value  of  the  semi-major  axis  'a'  of  the  adopted  ideal  reference  ellipsoid  obtained  as 
a=63781 36.60  iK).05  m,  for  both  the  Modified  Stokes'  and  the  LSC  technique,  as  a  by¬ 
product  to  this  study  agrees  well  with  the  recent  estimates  from  TOPEX/Poseidon  data  (Rapp, 
1993,  personal  communication). 

6.  Estimation  of  the  separation  between  various  regional  vertical  surfaces  from  the  realized  GVD 
was  done  on  the  assumption  that  the  regional  vertical  datums  were  established  to  the  required 
accuracy,  which  in  turn  can  help  in  detining  the  GVD  to  ±10  accuracy  level.  In  many  of  the 
regional  vertical  datums  only  a  single  station  has  been  used  assuming  that  it  represents  the  local 
height  system  in  entirety,  which  may  not  be  true.  Also  bias  caused  by  using  ^vity  anomalies 
wiA  certain  systematic  errors  (for  example,  gravity  anomalies  not  corrected  for  terrain)  should 
also  be  kept  in  mind  when  looking  at  the  estimated  values  of  separation  between  regional 
surfaces. 

7.  The  estimates  of  the  parameters  defining  the  GVD  computed  both  from  riKxiified  Stokes'  and 
LSC  techniques  are  generally  in  agreement,  except  for  some  numerical  differences.  These 
numerical  differences  are  mainly  due  to  the  differences  in  computed  gravimetric  height 
anomaly/undulation  between  the  two  techniques  which  varied  from  20  to  50  cm  at  certain 
stations.  Since  we  use  homogeneous  and  dense  gravity  data  around  the  computation  point,  one 
can  consider  that  for  practical  computation  of  height  anomaly/undulation,  modifi^  Stokes' 
technique  provides  a  better  solution.  As  also  expressed  by  Moritz  (1976,  p.  38),  even  if  the 
data  distribution  is  not  ideal,  it  is  better  to  use  collocation  for  interpolation,  in  order  to  densify 
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the  data  distribution,  and  subsequently  apply  Stokes’  method,  rather  than  using  LSC  technique 
for  a  straight  one  step  solutitm. 

8.  Finally,  tlM  results  reported  here  are  promising  enough  to  warrant  actual  implementatitm  of  a 
Global  Vertical  Datum  that  can  help  one  go  from  one  regional  vertical  datum  to  another 
knowing  the  single  bias  function  that  separates  the  two  datums. 
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CHAPTER? 


RECOMMENDATIONS  FOR  THE  ACQUISITION  AND  USE  OF  DATA  FOR  DEFINING 
AND  REALEING  A  GLOBAL  VERTICAL  DATUM 


In  Section  3.1,  the  different  kinds  of  data  that  are  essential  for  defining  and  realizing  a 
Global  Vertical  Datum  (GVD)  are  discussed  in  (tetail.  A  recommended  strategy  for  acquiring  and 
use  of  such  data  for  the  implementation  of  a  GVD  are  discussed  in  this  chapter.  The  strategy 
recommended  is  conceptually  simple  and  built  on  ongdng  International  activities. 

The  definition  and  realization  of  the  GVD  should  be  developed  around  the  lERS  terrestrial 
reference  frame  using  the  precise  positioning  data  given  by  the  SLR/VLBI  observations  all  over  the 
world.  There  are  about  171  lERS  sites  (refer  to  figure  7.1)  carrying  out  noobile  and  fixed 
SLR/VLBI  observations  distributed  in  over  42  countries  (Boucher  et  al.,  1993,  Table  1).  In 
addition  there  are  85  operational/proposed  IGS  stations  (refer  to  figure  7.2)  and  several  DORIS 
tracking  stations  (refer  to  figure  7.3),  adequately  complementing  the  SLR/VLBI  stations  forming  a 
dense  mesh  of  space  geodetic  stations  which  can  provide  the  required  precise  ellipsoidal  heights  for 
the  realization  of  G^.  Global  geopotential  tnodels  to  degtM  360  are  currently  available  and 
efforts  are  also  on  to  inqirove  the  accuracy  of  such  tnodels. 

As  the  problem  of  implementing  a  Global  Vertical  Datum,  acceptable  on  a  global  basis,  is 
of  longterm  character  requiring  close  cooperation  between  specialists  from  various  countries,  a 
special  study  group  similar  to  lAG  SSG  1.75  (Rrq^,  1987)  under  aegis  of  lAG  can  be  recreated  or 
the  objective  of  the  currently  operational  lAG  SSG  5.149  on  'Studies  on  Vertical  datums'  can  be 
enhanced  with  more  members  from  participating  countries.  This  SSG  in  additicm  to  identifying  a 
suitable  International  organization  to  coordmte  Ae  activities  of  the  various  scientific  organizations 
in  (Ufferent  countries,  >^1  also  lay  down  the  procedure  for  data  collection,  storage,  manipulation 
and  transfer  of  data  and  results  by  the  participating  countries. 

The  individual  participating  countries  should  be  made  responsible  for  the  accurate 
computaticMi  of  gravimetric  geoid  undulation  using  a  dense  grid  of  anotrialy  and  terrain  data  around 
the  computation  point  and  also  the  same  geopotential  rrxxiel  considered  to  be  accurate  and  used  by 
other  participating  countries.  They  will  also  provide  the  International  organizations  comdinating 
this  collective  endeavOT  with  the  precise  height  of  the  space  geodetic  stations,  referring  to  the  same 
point  to  which  ellipsoidal  height  is  referred,  following  the  1st  order  leveling  procedures  from  an 
adjusted  Bench  Marie  in  the  regional  vertical  refermce  frame. 

Finally,  by  building  on  existing  structures,  adding  to  the  existing  networks,  and  adapting 
the  existing  infrastructure  to  the  expanded  needs  associated  with  the  Global  Vertical  Datum 
concept,  we  can  achieve  a  gradual  and  flexible  implementation.  For  this  we  need  to  begin  the 
operation  of  realizing  a  G^  immediately,  using  existing  data  sources.  The  development  of 
r^uired  mathematic^  models  and  the  first  iteration  attempt  in  realizing  a  GVD  discussed  in  this 
dissertation  report  are  the  first  few  steps  in  that  direction. 
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Figure  7. 1  Map  of  Geodetic  Sites  -  Global  (NASA  Tech.  Memo  4482) 
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Figure  7.2  GPS  Tracking  Network  of  the  International  GPS  Service  for  Geodynamics  - 
Operational  and  Proposed  Stations 
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Figure  7.3  DORIS  Orbitography  Beacon  Network  (Dec.  1990) 
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APPENDIX  A 


COMPUTATION  OF  GRAVIMETRIC  HEIGHT  ANOMALYAJNDULATION  USING  JGM-2 
(AUGMENTED)  GEOPOTENTIAL  MODEL  (TO  DEGREE  360) 

The  OSU91A  geopotential  model  which  was  used  as  the  reference  model  in  the 
numerical  investigadons  reported  in  Chapter  4  was  replaced  by  the  newly  made  available  JGM-2 
(NASA  model)  (Nerem  et  al.,  1993)  fr^  degree  2  to  70  and  augmented  by  OSU91A  model 
coefficients  from  degree  71  to  360,  and  the  computation  of  gravimetric  height  anomaly  Cl  and 
undulation  Nj  fiom  equations  (3.71)  and  (3.72)  was  repeated  for  all  space  geodetic  stadcms  used 
in  the  GVD  defmition  study.  Table  Al  shows  the  newly  computed  vtdues  of  remote  zone 
contribution  to  gravimetric  height  anomalyAinduladon  computation  using  JGM-2  (augmented 
model)  in  column  (2),  with  ciq>  contribudon  in  column  (1)  and  other  conecdons  in  column  (3) 
reproduced  from  Tables  (4.7)  to  (4.9b)  along  widi  the  odier  components  reproduced  frcmi  Tables 
(4.7),  (4.8),  (4.9a)  and  (4.9b).  Total  height  anomaly/unduladon  as  well  as  the  Y  values  computed 
using  OSU91 A  model  and  JGM-2  (augmented)  inodel  are  also  shown  for  better  comparison  of 
the  variadon  in  C|/N|  computadon  using  different  global  geopotendal  model  coefficients. 

From  Table  Al  we  see  that  the  height  anomalyAmduladons  computed  using  the  OSU91A 
reference  model  differ  from  the  one  computed  using  JGM-2  (augmented)  model  by  about  2  to  17 
cm  at  the  space  geodetic  stations  in  the  United  States,  4  to  2S  cm  at  the  European  and 
Scandinavian  stations,  and  about  S  to  9  cm  at  the  Australian  stations.  Since  in  the  JGM-2 
(augmented)  model  die  lower  degree  harmonics  have  been  computed  more  accurately  than  the 
OSU91A  model,  the  difference  noted  in  the  ctmqiuted  undulations  can  be  considered  to  be  due  to 
long  wavelength  error  in  the  OSU91A  geopmential  model  used  as  the  reference  model  in  the 
numerical  stupes.  Ilie  differences  between  the  height  anomaly^ndulati(Mls  cmi^uted  using  the 
OSU91A  reference  model  and  the  JGM-2  (augmented)  model  are  consistent  with  the  accuracy 
estimate  given  in  Table  4.19  for  the  Modifi^  Stdces’  technique.  The  estimated  parameters 
deHning  die  Global  Vertical  Datum  computed  using  the  Y-values  from  JGM-2  (augnwnted) 
model  are  shown  in  Figure  Al.  Comparison  of  Figures  4.S  and  Al  show  that  the  estimated 
parameters  agree  to  an  extent  of  ±10  cm  in  all  cases  except  for  ODN  datum  where  the 
disagreement  is  about  19  cm. 


NAVD88 


ODN 


Figure  A1  Regicmal  Veitical  Datum  Separatitm  From  Defined  Global  Vertical  Datum  (computed 
using  JGM-2  (augmented)  Model  and  the  Modified  Stokes'  Technique).  Units  are  in 


Table  A1 :  Comparison  of  Y-values  Computed  Using  Different  Global  Geopotential  Models. 
Units  are  in  meters. 
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